Change in Backward Euler Feap 8.4 and 8.6

[user2705]

Dear,

Did something change in implementation of backward Euler in Feap 8.6 from Feap 8.4?

Best regards and thanks in advance

[FEAP_Admin]

Yes.  Part of the code has been re-written with respect to the rate vectors and the argument list to program/dyna02.f has changed.
Are you seeing an error?

Here is a diff of the two routines: > is ver84  < is ver86

Code: [Select]

diff dyna02.f ../../ver84/program/
1,2c1,2
< !$Id:$
<       subroutine dyna02(du,urate,nneq,ndf,ndp,ndo,npart,isw)
---
> c$Id:$
>       subroutine dyna02(du,urate,nneq,ndf,ndp,ndo,isw)
4c4
< !      * * F E A P * * A Finite Element Analysis Program
---
> c      * * F E A P * * A Finite Element Analysis Program
6,7c6,7
< !....  Copyright (c) 1984-2024: Regents of the University of California
< !                               All rights reserved
---
> c....  Copyright (c) 1984-2014: Regents of the University of California
> c                               All rights reserved
9,41c9,38
< !-----[--.----+----.----+----.-----------------------------------------]
< !     Modification log                                Date (dd/mm/year)
< !       Original version                                    01/11/2006
< !       1. Move backup storage to 'nrt' based index         20/10/2010
< !       2. Add 'npart' to argument delete 'part0.h'         23/11/2019
< !       3. Modify update for nod(*) = 2                     02/02/2020
< !-----[--.----+----.----+----.-----------------------------------------]
< !      Purpose: Static, 1st and 2nd order ODE integration by backward
< !               Euler method.
<
< !      Notes: 1. Values of velocity and acceleration are returned as
< !                zero for values of ndo(i) .gt. order specified.
<
< !      Inputs:
< !         du(*)             Increment to displacement
< !         urate(nneq,*)     Rate vectors - fixed by ALGO
< !         nneq              numnp * ndf
< !         ndf               Number of DOF/node
< !         ndp(*)            Partition dof's
< !         ndo(*)            Order dof's
< !         npart             Active partition number
< !         isw               Control switch
< !                            1  STARTING update: begining of time step
< !                            2  UPDATE at an iteration within time step
< !                            3  BACK solution to begining of time step
<
< !      Outputs:
< !         urate(nneq,nn)    Rate vectors:
< !                            1  Velocity at t_n+1     (ndo(i) .ge. 1
< !                            2  Acceleration at t_n+1 (ndo(i) .ge. 2
< !                            6  Velocity at t_n       (ndo(i) .ge. 1
< !                            7  Acceleration at t_n   (ndo(i) .ge. 2
< !-----[--.----+----.----+----.-----------------------------------------]
---
> c-----[--.----+----.----+----.-----------------------------------------]
> c     Modification log                                Date (dd/mm/year)
> c       Original version                                    01/11/2006
> c       1. Move backup storage to 'nrt' based index         20/10/2010
> c-----[--.----+----.----+----.-----------------------------------------]
> c      Purpose: Static, 1st and 2nd order ODE integration by backward
> c               Euler method.
>
> c      Notes: 1. Values of velocity and acceleration are returned as
> c                zero for values of ndo(i) .gt. order specified.
>
> c      Inputs:
> c         du(*)             Increment to displacement
> c         urate(nneq,*)     Rate vectors - fixed by ALGO
> c         nneq              numnp * ndf
> c         ndf               Number of DOF/node
> c         ndp(*)            Partition dof's
> c         ndo(*)            Order dof's
> c         isw               Control switch
> c                            1  STARTING update: begining of time step
> c                            2  UPDATE at an iteration within time step
> c                            3  BACK solution to begining of time step
>
> c      Outputs:
> c         urate(nneq,nn)    Rate vectors:
> c                            1  Velocity at t_n+1     (ndo(i) .ge. 1
> c                            2  Acceleration at t_n+1 (ndo(i) .ge. 2
> c                            6  Velocity at t_n       (ndo(i) .ge. 1
> c                            7  Acceleration at t_n   (ndo(i) .ge. 2
> c-----[--.----+----.----+----.-----------------------------------------]
44a42
>       include  'part0.h'
47,48c45,46
<       integer       :: i, n, nneq,ndf,isw, ndp(*),ndo(*),npart
<       real (kind=8) :: du(*),urate(nneq,*)
---
>       integer   i, n, nneq,ndf,isw, ndp(*),ndo(*)
>       real*8    du(*),urate(nneq,*)
52c50,51
< !     Backward Euler: Initialize at start of step
---
> c     Backward Euler: Initialize at start of step
>
63,64c62,63
<                 urate(n,2) = -c1*urate(n,1)
<                 urate(n,1) =  0.0d0
---
>                 urate(n,1) = -   urate(n,1)
>                 urate(n,2) =     urate(n,1)
74c73,74
< !     Backward Euler: Updates in iterations
---
> c     Backward Euler: Updates in iterations
>
82c82
<                 urate(n,2) = urate(n,2) + c1*urate(n,1)
---
>                 urate(n,2) = urate(n,1)
93c93,94
< !     Backward Euler: Back up to start of step
---
> c     Backward Euler: Back up to start of step
>
107c108
<       end subroutine dyna02
---
>       end



[user2705]

Dear Feap Admin,

Thanks for the reply.

I had an element in Feap 8.4. that had to have small mass in the model in order to work. Because my problem was first order, I was solving it with TRANS BACK, since in FEAP 8.4. this scheme was accommodated to solve second order problems (I had a previous post about this topic on forum).

In Feap 8.6. this element cannot be started with TRANS BACK. For me it is ok to use TRANS NEWM, I was just curious to see whether the implementation of backward Euler has changed, particularly its implementation regarding second order problems.

Best regards and thanks in advance

[Prof. R.L. Taylor]

There is an error in version 8.6 for backward Euler applied to problems with accelerations (d^2 u/dt^2) terms.

Correctioins are:

1. In ./program/dsetci.f:  set gtan(3) = c1*c1

2. In ./program/dyna02.f:  in the update after isw.eq.2 set: urate(n,2) = urate(n,2) + c21*c1*du(n)

be sure to check that a linear problem converges in the second iteration.

[user2705]

Dear professor,

Thank you very much for the reply.

Did you mean here urate(n,2) = urate(n,2) + c21*c1*du(n) instead of c21--c1?

I also have one more question regarding generalized mid point scheme for first order problems. My model works fine with TRANS BACK and TRANS NEWMARK, but it does not converge with default TRANS GEN1 (with alpha=1 is ok). What could be the issue?

Best regards and thanks in advance

[Prof. R.L. Taylor]

Sorry should be c1*c1

[user2705]

Dear professor,

Thank you very much.

Regarding my second question, I need to intervene in my user element.

Best regards

[Prof. R.L. Taylor]

I am not sure GEN1 is coded for transients.  Look at the manual.  You would need your stiffness tangent multiplied by ctan(1) for it to work.  Look at dyna07.f to see how it works.  Does HHT work?

[user2705]

Dear professor,

It was my mistake. I did not realize immediately that I will have to modify my user element for trans,gen1,0.5 to work for the first order problem I am running, and energy conserving scheme for another second order problem (with this element in which I need to have small mass).

Best regards and thank you very much for your kind help