CalculiX  2.13
A Free Software Three-Dimensional Structural Finite Element Program
resultsini.f File Reference

Go to the source code of this file.

Functions/Subroutines

subroutine resultsini (nk, v, ithermal, filab, iperturb, f, fn, nactdof, iout, qa, vold, b, nodeboun, ndirboun, xboun, nboun, ipompc, nodempc, coefmpc, labmpc, nmpc, nmethod, cam, neq, veold, accold, bet, gam, dtime, mi, vini, nprint, prlab, intpointvarm, calcul_fn, calcul_f, calcul_qa, calcul_cauchy, nener, ikin, intpointvart, xforc, nforc)
 

Function/Subroutine Documentation

◆ resultsini()

subroutine resultsini ( integer  nk,
real*8, dimension(0:mi(2),*)  v,
integer, dimension(2)  ithermal,
character*87, dimension(*)  filab,
integer, dimension(*)  iperturb,
real*8, dimension(*)  f,
real*8, dimension(0:mi(2),*)  fn,
integer, dimension(0:mi(2),*)  nactdof,
integer  iout,
real*8, dimension(*)  qa,
real*8, dimension(0:mi(2),*)  vold,
real*8, dimension(*)  b,
integer, dimension(*)  nodeboun,
integer, dimension(*)  ndirboun,
real*8, dimension(*)  xboun,
integer  nboun,
integer, dimension(*)  ipompc,
integer, dimension(3,*)  nodempc,
real*8, dimension(*)  coefmpc,
character*20, dimension(*)  labmpc,
integer  nmpc,
integer  nmethod,
real*8, dimension(5)  cam,
integer  neq,
real*8, dimension(0:mi(2),*)  veold,
real*8, dimension(0:mi(2),*)  accold,
real*8  bet,
real*8  gam,
real*8  dtime,
integer, dimension(*)  mi,
real*8, dimension(0:mi(2),*)  vini,
integer  nprint,
character*6, dimension(*)  prlab,
integer  intpointvarm,
integer  calcul_fn,
integer  calcul_f,
integer  calcul_qa,
integer  calcul_cauchy,
integer  nener,
integer  ikin,
integer  intpointvart,
real*8, dimension(*)  xforc,
integer  nforc 
)
25 !
26 ! initialization
27 !
28 ! 1. storing the calculated primary variables nodewise
29 ! 2. inserting the boundary conditions nodewise (SPC's and MPC's)
30 ! 3. determining which derived variables (strains, stresses,
31 ! internal forces...) have to be calculated
32 !
33  implicit none
34 !
35  character*6 prlab(*)
36  character*20 labmpc(*)
37  character*87 filab(*)
38 !
39  integer mi(*),nactdof(0:mi(2),*),nodeboun(*),ndirboun(*),
40  & ipompc(*),nodempc(3,*),mt,nk,ithermal(2),i,j,
41  & nener,iperturb(*),iout,nboun,nmpc,nmethod,ist,ndir,node,index,
42  & neq,incrementalmpc,nprint,ikin,calcul_fn,nforc,
43  & calcul_f,calcul_cauchy,calcul_qa,intpointvarm,intpointvart,
44  & irefnode,irotnode,iexpnode,irefnodeprev
45 !
46  real*8 v(0:mi(2),*),vini(0:mi(2),*),f(*),fn(0:mi(2),*),
47  & cam(5),vold(0:mi(2),*),b(*),xboun(*),coefmpc(*),
48  & veold(0:mi(2),*),accold(0:mi(2),*),xforc(*),
49  & qa(*),bet,gam,dtime,scal1,scal2,bnac,
50  & fixed_disp
51 !
52  mt=mi(2)+1
53 !
54  if((iout.ne.2).and.(iout.gt.-1)) then
55 !
56  if((nmethod.ne.4).or.(iperturb(1).le.1)) then
57  if(ithermal(1).ne.2) then
58  do i=1,nk
59  do j=1,mi(2)
60  if(nactdof(j,i).gt.0) then
61  bnac=b(nactdof(j,i))
62  else
63  cycle
64  endif
65  v(j,i)=v(j,i)+bnac
66  if((iperturb(1).ne.0).and.(abs(nmethod).eq.1)) then
67  if(dabs(bnac).gt.cam(1)) then
68  cam(1)=dabs(bnac)
69  cam(4)=nactdof(j,i)-0.5d0
70  endif
71  endif
72  enddo
73  enddo
74  endif
75  if(ithermal(1).gt.1) then
76  do i=1,nk
77  if(nactdof(0,i).gt.0) then
78  bnac=b(nactdof(0,i))
79  else
80  cycle
81  endif
82  v(0,i)=v(0,i)+bnac
83  if((iperturb(1).ne.0).and.(abs(nmethod).eq.1)) then
84  if(dabs(bnac).gt.cam(2)) then
85  cam(2)=dabs(bnac)
86  cam(5)=nactdof(0,i)-0.5d0
87  endif
88  endif
89  enddo
90  endif
91 !
92  else
93 !
94 ! direct integration dynamic step
95 ! b contains the acceleration increment
96 !
97  if(ithermal(1).ne.2) then
98  scal1=bet*dtime*dtime
99  scal2=gam*dtime
100  do i=1,nk
101  do j=1,mi(2)
102  if(nactdof(j,i).gt.0) then
103  bnac=b(nactdof(j,i))
104  else
105  cycle
106  endif
107  v(j,i)=v(j,i)+scal1*bnac
108  if(dabs(scal1*bnac).gt.cam(1)) then
109  cam(1)=dabs(scal1*bnac)
110  cam(4)=nactdof(j,i)-0.5d0
111  endif
112  veold(j,i)=veold(j,i)+scal2*bnac
113  accold(j,i)=accold(j,i)+bnac
114  enddo
115  enddo
116  endif
117  if(ithermal(1).gt.1) then
118  do i=1,nk
119 !
120 ! no extrapolation is performed for the solution
121 ! of parabolic differential equations (e.g. heat
122 ! equation)
123 !
124  veold(0,i)=0.d0
125  if(nactdof(0,i).gt.0) then
126  bnac=b(nactdof(0,i))
127  else
128  cycle
129  endif
130  v(0,i)=v(0,i)+bnac
131  if(dabs(bnac).gt.cam(2)) then
132  cam(2)=dabs(bnac)
133  cam(5)=nactdof(0,i)-0.5d0
134  endif
135  cam(3)=max(cam(3),dabs(v(0,i)-vini(0,i)))
136  enddo
137  endif
138  endif
139 !
140  endif
141 !
142 ! initialization
143 !
144  calcul_fn=0
145  calcul_f=0
146  calcul_qa=0
147  calcul_cauchy=0
148 !
149 ! determining which quantities have to be calculated
150 !
151  if((iperturb(1).ge.2).or.((iperturb(1).le.0).and.(iout.lt.0)))
152  & then
153  if((iout.lt.1).and.(iout.gt.-2)) then
154  calcul_fn=1
155  calcul_f=1
156  calcul_qa=1
157  elseif((iout.ne.-2).and.(iperturb(2).eq.1)) then
158  calcul_cauchy=1
159  endif
160  endif
161 !
162  if(iout.gt.0) then
163  if((filab(5)(1:4).eq.'RF ').or.
164  & (filab(10)(1:4).eq.'RFL ')) then
165  calcul_fn=1
166  else
167  do i=1,nprint
168  if((prlab(i)(1:4).eq.'RF ').or.
169  & (prlab(i)(1:4).eq.'RFL ')) then
170  calcul_fn=1
171  exit
172  endif
173  enddo
174  endif
175  endif
176 !
177 ! initializing fn
178 !
179  if(calcul_fn.eq.1) then
180  do i=1,nk
181  do j=0,mi(2)
182  fn(j,i)=0.d0
183  enddo
184  enddo
185  endif
186 !
187 ! initializing f
188 !
189  if(calcul_f.eq.1) then
190  do i=1,neq
191  f(i)=0.d0
192  enddo
193  endif
194 !
195 ! SPC's and MPC's have to be taken into account for
196 ! iout=0,1 and -1
197 !
198  if(abs(iout).lt.2) then
199 !
200 ! inserting the boundary conditions
201 !
202  do i=1,nboun
203  if(ndirboun(i).gt.mi(2)) cycle
204  fixed_disp=xboun(i)
205 !
206 ! a discontinuity in the displacements in an initial
207 ! acceleration step (recognized by the "special" time
208 ! increment) should not lead to a change in
209 ! the acceleration or velocity; actually, such a
210 ! discontinuity is not allowed since it leads to
211 ! infinite accelerations
212 !
213 c if((nmethod.eq.4).and.(iperturb(1).gt.1).and.
214 c & (nint(dtime*1.d28).ne.123571113)) then
215  if((nmethod.eq.4).and.(iperturb(1).gt.1)) then
216  ndir=ndirboun(i)
217  node=nodeboun(i)
218  if(ndir.gt.0) then
219 !
220 ! bnac is the change of acceleration
221 !
222  bnac=(xboun(i)-v(ndir,node))/
223  & (bet*dtime*dtime)
224  if(nint(dtime*1.d28).ne.123571113) then
225  veold(ndir,node)=veold(ndir,node)+
226  & gam*dtime*bnac
227  endif
228  accold(ndir,node)=accold(ndir,node)+bnac
229  endif
230  endif
231  v(ndirboun(i),nodeboun(i))=fixed_disp
232  enddo
233 !
234 ! inserting the mpc information
235 ! the parameter incrementalmpc indicates whether the
236 ! incremental displacements enter the mpc or the total
237 ! displacements (incrementalmpc=0)
238 !
239 c
240 c to be checked: should replace the lines underneath do i=1,nmpc
241 c
242 c incrementalmpc=iperturb(2)
243 !
244  do i=1,nmpc
245  if((labmpc(i)(1:20).eq.' ').or.
246  & (labmpc(i)(1:6).eq.'CYCLIC').or.
247  & (labmpc(i)(1:9).eq.'SUBCYCLIC')) then
248  incrementalmpc=0
249  else
250  if((nmethod.eq.2).or.(nmethod.eq.3).or.
251  & ((iperturb(1).eq.0).and.(abs(nmethod).eq.1)))
252  & then
253  incrementalmpc=0
254  else
255  incrementalmpc=1
256  endif
257  endif
258  ist=ipompc(i)
259  node=nodempc(1,ist)
260  ndir=nodempc(2,ist)
261  if(ndir.eq.0) then
262  if(ithermal(1).lt.2) cycle
263  elseif(ndir.gt.mi(2)) then
264  cycle
265  else
266  if(ithermal(1).eq.2) cycle
267  endif
268  index=nodempc(3,ist)
269  fixed_disp=0.d0
270  if(index.ne.0) then
271  do
272  if(incrementalmpc.eq.0) then
273  fixed_disp=fixed_disp-coefmpc(index)*
274  & v(nodempc(2,index),nodempc(1,index))
275  else
276  fixed_disp=fixed_disp-coefmpc(index)*
277  & (v(nodempc(2,index),nodempc(1,index))-
278  & vold(nodempc(2,index),nodempc(1,index)))
279  endif
280  index=nodempc(3,index)
281  if(index.eq.0) exit
282  enddo
283  endif
284  fixed_disp=fixed_disp/coefmpc(ist)
285  if(incrementalmpc.eq.1) then
286  fixed_disp=fixed_disp+vold(ndir,node)
287  endif
288 !
289 ! a discontinuity in the displacements in an initial
290 ! acceleration step (recognized by the "special" time
291 ! increment) should not lead to a change in
292 ! the acceleration or velocity; actually, such a
293 ! discontinuity is not allowed since it leads to
294 ! infinite accelerations
295 !
296 c if((nmethod.eq.4).and.(iperturb(1).gt.1).and.
297 c & (nint(dtime*1.d28).ne.123571113)) then
298  if((nmethod.eq.4).and.(iperturb(1).gt.1)) then
299  if(ndir.gt.0) then
300 !
301 ! bnac is the change of acceleration
302 !
303  bnac=(fixed_disp-v(ndir,node))/
304  & (bet*dtime*dtime)
305  if(nint(dtime*1.d28).ne.123571113) then
306  veold(ndir,node)=veold(ndir,node)+
307  & gam*dtime*bnac
308  endif
309  accold(ndir,node)=accold(ndir,node)+bnac
310  endif
311  endif
312  v(ndir,node)=fixed_disp
313  enddo
314  endif
315 !
316 ! storing the knot information in the .dat-file
317 !
318  if(ithermal(1).ne.2) then
319  irefnodeprev=0
320  do i=1,nmpc
321  if(iout.gt.0) then
322  if(labmpc(i)(1:4).eq.'KNOT') then
323  irefnode=nodempc(1,nodempc(3,ipompc(i)))
324  if(irefnode.ne.irefnodeprev) then
325  irefnodeprev=irefnode
326  iexpnode=nodempc(1,nodempc(3,nodempc(3,ipompc(i))))
327  if(labmpc(i)(5:5).ne.'2') then
328  irotnode=nodempc(1,nodempc(3,nodempc(3,
329  & nodempc(3,ipompc(i)))))
330  else
331  irotnode=nodempc(1,nodempc(3,nodempc(3,
332  & nodempc(3,nodempc(3,nodempc(3,ipompc(i)))))))
333  endif
334 c write(5,*)
335 c write(5,'(a5)') labmpc(i)(1:5)
336 c write(5,'("tra",i10,3(1x,e11.4))')
337 c & irefnode,(v(j,irefnode),j=1,3)
338 c write(5,'("rot",i10,3(1x,e11.4))')
339 c & irotnode,(v(j,irotnode),j=1,3)
340 c if(labmpc(i)(5:5).eq.'2') then
341 c write(5,'("exp",i10,3(1x,e11.4))')
342 c & iexpnode,(v(j,iexpnode),j=1,3)
343 c else
344 c write(5,'("exp",i10,3(1x,e11.4))')
345 c & iexpnode,v(1,iexpnode)
346 c endif
347  endif
348  endif
349  endif
350  enddo
351  endif
352 !
353 ! check whether there are any strain output requests
354 !
355  ikin=0
356 !
357 ! *DYNAMIC
358 !
359  if((nmethod.eq.4).and.(iperturb(1).gt.1).and.
360  & (ithermal(1).le.1)) then
361  ikin=1
362  endif
363 !
364 ! dat-output
365 !
366  do i=1,nprint
367  if((prlab(i)(1:4).eq.'ENER').or.(prlab(i)(1:4).eq.'ELSE').or.
368  & (prlab(i)(1:4).eq.'CELS')) then
369  elseif(prlab(i)(1:4).eq.'ELKE') then
370  ikin=1
371  endif
372  enddo
373 !
374  qa(1)=0.d0
375  qa(2)=0.d0
376 !
377 ! check whether integration point variables are needed in
378 ! modal dynamics and steady state dynamics calculations
379 !
380  intpointvarm=1
381  intpointvart=1
382 !
383  if((nmethod.ge.4).and.(nmethod.le.5).and.(iperturb(1).lt.2)) then
384  intpointvarm=0
385  if((filab(3)(1:4).eq.'S ').or.
386  & (filab(4)(1:4).eq.'E ').or.
387  & (filab(5)(1:4).eq.'RF ').or.
388  & (filab(6)(1:4).eq.'PEEQ').or.
389  & (filab(7)(1:4).eq.'ENER').or.
390  & (filab(8)(1:4).eq.'SDV ').or.
391  & (filab(13)(1:4).eq.'ZZS ').or.
392  & (filab(13)(1:4).eq.'ERR ').or.
393  & (filab(18)(1:4).eq.'PHS ').or.
394  & (filab(20)(1:4).eq.'MAXS').or.
395  & (filab(26)(1:4).eq.'CONT').or.
396  & (filab(27)(1:4).eq.'CELS')) intpointvarm=1
397  do i=1,nprint
398  if((prlab(i)(1:4).eq.'S ').or.
399  & (prlab(i)(1:4).eq.'E ').or.
400  & (prlab(i)(1:4).eq.'PEEQ').or.
401  & (prlab(i)(1:4).eq.'ENER').or.
402  & (prlab(i)(1:4).eq.'ELKE').or.
403  & (prlab(i)(1:4).eq.'CDIS').or.
404  & (prlab(i)(1:4).eq.'CSTR').or.
405  & (prlab(i)(1:4).eq.'CELS').or.
406  & (prlab(i)(1:4).eq.'SDV ').or.
407  & (prlab(i)(1:4).eq.'RF ')) then
408  intpointvarm=1
409  exit
410  endif
411  enddo
412 !
413  intpointvart=0
414  if((filab(9)(1:4).eq.'HFL ').or.
415  & (filab(10)(1:4).eq.'RFL ')) intpointvart=1
416  do i=1,nprint
417  if((prlab(i)(1:4).eq.'HFL ').or.
418  & (prlab(i)(1:4).eq.'RFL ')) intpointvart=1
419  enddo
420 !
421 ! if internal forces are requested integration point
422 ! values have to be calculated
423 !
424  if(calcul_fn.eq.1) then
425  intpointvarm=1
426  intpointvart=1
427  endif
428  endif
429 !
430  return
#define max(a, b)
Definition: cascade.c:32
Hosted by OpenAircraft.com, (Michigan UAV, LLC)