CalculiX  2.8
A Free Software Three-Dimensional Structural Finite Element Program
 All Classes Files Functions Variables Macros
checkconvergence.c File Reference
#include <stdio.h>
#include <math.h>
#include <stdlib.h>
#include "CalculiX.h"
Include dependency graph for checkconvergence.c:

Go to the source code of this file.

Functions

void checkconvergence (double *co, ITG *nk, ITG *kon, ITG *ipkon, char *lakon, ITG *ne, double *stn, ITG *nmethod, ITG *kode, char *filab, double *een, double *t1act, double *time, double *epn, ITG *ielmat, char *matname, double *enern, double *xstaten, ITG *nstate_, ITG *istep, ITG *iinc, ITG *iperturb, double *ener, ITG *mi, char *output, ITG *ithermal, double *qfn, ITG *mode, ITG *noddiam, double *trab, ITG *inotr, ITG *ntrans, double *orab, ITG *ielorien, ITG *norien, char *description, double *sti, ITG *icutb, ITG *iit, double *dtime, double *qa, double *vold, double *qam, double *ram1, double *ram2, double *ram, double *cam, double *uam, ITG *ntg, double *ttime, ITG *icntrl, double *theta, double *dtheta, double *veold, double *vini, ITG *idrct, double *tper, ITG *istab, double *tmax, ITG *nactdof, double *b, double *tmin, double *ctrl, double *amta, ITG *namta, ITG *itpamp, ITG *inext, double *dthetaref, ITG *itp, ITG *jprint, ITG *jout, ITG *uncoupled, double *t1, ITG *iitterm, ITG *nelemload, ITG *nload, ITG *nodeboun, ITG *nboun, ITG *itg, ITG *ndirboun, double *deltmx, ITG *iflagact, char *set, ITG *nset, ITG *istartset, ITG *iendset, ITG *ialset, double *emn, double *thicke, char *jobnamec, ITG *mortar, ITG *nmat, ITG *ielprop, double *prop)
 

Function Documentation

void checkconvergence ( double *  co,
ITG nk,
ITG kon,
ITG ipkon,
char *  lakon,
ITG ne,
double *  stn,
ITG nmethod,
ITG kode,
char *  filab,
double *  een,
double *  t1act,
double *  time,
double *  epn,
ITG ielmat,
char *  matname,
double *  enern,
double *  xstaten,
ITG nstate_,
ITG istep,
ITG iinc,
ITG iperturb,
double *  ener,
ITG mi,
char *  output,
ITG ithermal,
double *  qfn,
ITG mode,
ITG noddiam,
double *  trab,
ITG inotr,
ITG ntrans,
double *  orab,
ITG ielorien,
ITG norien,
char *  description,
double *  sti,
ITG icutb,
ITG iit,
double *  dtime,
double *  qa,
double *  vold,
double *  qam,
double *  ram1,
double *  ram2,
double *  ram,
double *  cam,
double *  uam,
ITG ntg,
double *  ttime,
ITG icntrl,
double *  theta,
double *  dtheta,
double *  veold,
double *  vini,
ITG idrct,
double *  tper,
ITG istab,
double *  tmax,
ITG nactdof,
double *  b,
double *  tmin,
double *  ctrl,
double *  amta,
ITG namta,
ITG itpamp,
ITG inext,
double *  dthetaref,
ITG itp,
ITG jprint,
ITG jout,
ITG uncoupled,
double *  t1,
ITG iitterm,
ITG nelemload,
ITG nload,
ITG nodeboun,
ITG nboun,
ITG itg,
ITG ndirboun,
double *  deltmx,
ITG iflagact,
char *  set,
ITG nset,
ITG istartset,
ITG iendset,
ITG ialset,
double *  emn,
double *  thicke,
char *  jobnamec,
ITG mortar,
ITG nmat,
ITG ielprop,
double *  prop 
)

Definition at line 33 of file checkconvergence.c.

53  {
54 
55  ITG i0,ir,ip,ic,il,ig,ia,iest,iest1=0,iest2=0,iconvergence,idivergence,
56  ngraph=1,k,*ipneigh=NULL,*neigh=NULL,*inum=NULL,id,istart,iend,inew,
57  i,j,mt=mi[1]+1,iexceed;
58 
59  double df,dc,db,dd,ran,can,rap,ea,cae,ral,da,*vr=NULL,*vi=NULL,*stnr=NULL,
60  *stni=NULL,*vmax=NULL,*stnmax=NULL,*cs=NULL,c1[2],c2[2],reftime,
61  *fn=NULL,*eenmax=NULL,*fnr=NULL,*fni=NULL,*qfx=NULL,*cdn=NULL,
62  *cdnr=NULL,*cdni=NULL;
63 
64  /* next lines are active if the number of contact elements was
65  changed in the present increment */
66 
67  if ((*iflagact==1)&&(*mortar==0)){
68  if(ctrl[0]<*iit+4)ctrl[0]=*iit+4;
69  if(ctrl[1]<*iit+8)ctrl[1]=*iit+8;
70  ctrl[3]+=1;
71  }
72 
73  i0=ctrl[0];ir=ctrl[1];ip=ctrl[2];ic=ctrl[3];il=ctrl[4];ig=ctrl[5];
74  ia=ctrl[7];df=ctrl[10];dc=ctrl[11];db=ctrl[12];da=ctrl[13];dd=ctrl[16];
75  ran=ctrl[18];can=ctrl[19];rap=ctrl[22];
76  ea=ctrl[23];cae=ctrl[24];ral=ctrl[25];
77 
78  /* check for forced divergence (due to divergence of a user material
79  routine */
80 
81  if(qa[2]>0.){idivergence=1;}else{idivergence=0;}
82 
83  if(*ithermal!=2){
84  if(qa[0]>ea*qam[0]){
85  if(*iit<=ip){c1[0]=ran;}
86  else{c1[0]=rap;}
87  c2[0]=can;
88  }
89  else{
90  c1[0]=ea;
91  c2[0]=cae;
92  }
93  if(ram1[0]<ram2[0]){ram2[0]=ram1[0];}
94  }
95  if(*ithermal>1){
96  if(qa[1]>ea*qam[1]){
97  if(*iit<=ip){c1[1]=ran;}
98  else{c1[1]=rap;}
99  c2[1]=can;
100  }
101  else{
102  c1[1]=ea;
103  c2[1]=cae;
104  }
105  if(ram1[1]<ram2[1]){ram2[1]=ram1[1];}
106  }
107 
108  iconvergence=0;
109 
110  /* mechanical */
111 
112  if(*ithermal<2){
113  if((*iit>1)&&(ram[0]<=c1[0]*qam[0])&&(*iflagact==0)&&
114 // if((*iit>1)&&(ram[0]<=c1[0]*qam[0])&&
115  ((cam[0]<=c2[0]*uam[0])||
116  (((ram[0]*cam[0]<c2[0]*uam[0]*ram2[0])||(ram[0]<=ral*qam[0])||
117  (qa[0]<=ea*qam[0]))&&(*ntg==0))||
118  (cam[0]<1.e-8))) iconvergence=1;
119  }
120 
121  /* thermal */
122 
123  if(*ithermal==2){
124  if((ram[1]<=c1[1]*qam[1])&&
125  (cam[2]<*deltmx)&&
126  ((cam[1]<=c2[1]*uam[1])||
127  (((ram[1]*cam[1]<c2[1]*uam[1]*ram2[1])||(ram[1]<=ral*qam[1])||
128  (qa[1]<=ea*qam[1]))&&(*ntg==0))||
129  (cam[1]<1.e-8)))iconvergence=1;
130  }
131 
132  /* thermomechanical */
133 
134  if(*ithermal==3){
135 // if(((ram[0]<=c1[0]*qam[0])&&
136  if(((*iit>1)&&(ram[0]<=c1[0]*qam[0])&&
137  ((cam[0]<=c2[0]*uam[0])||
138  (((ram[0]*cam[0]<c2[0]*uam[0]*ram2[0])||(ram[0]<=ral*qam[0])||
139  (qa[0]<=ea*qam[0]))&&(*ntg==0))||
140  (cam[0]<1.e-8)))&&
141  ((ram[1]<=c1[1]*qam[1])&&
142  (cam[2]<*deltmx)&&
143  ((cam[1]<=c2[1]*uam[1])||
144  (((ram[1]*cam[1]<c2[1]*uam[1]*ram2[1])||(ram[1]<=ral*qam[1])||
145  (qa[1]<=ea*qam[1]))&&(*ntg==0))||
146  (cam[1]<1.e-8))))iconvergence=1;
147  }
148 
149  /* reset iflagact */
150 
151  *iflagact=0;
152 
153  /* increment convergence reached */
154 
155  if((iconvergence==1)&&(idivergence==0)){
156 // *ttime=*ttime+*dtime;
157 
158  /* cutting the insignificant digits from ttime */
159 
160 // *ttime=*ttime+1.;
161 // *ttime=*ttime-1.;
162  FORTRAN(writesummary,(istep,iinc,icutb,iit,ttime,time,dtime));
163  if(*uncoupled){
164  if(*ithermal==2){
165  *iitterm=*iit;
166  *ithermal=1;
167  for(k=0;k<*nk;++k){t1[k]=vold[mt*k];}
168 // *ttime=*ttime-*dtime;
169  *iit=1;
170  (ctrl[0])*=4;
171  printf(" thermal convergence\n\n");
172  return;
173  }else{
174  *ithermal=3;
175  *iit=*iitterm;
176  (ctrl[0])/=4;
177  }
178  }
179 
180  *icntrl=1;
181  *icutb=0;
182  *theta=*theta+*dtheta;
183 
184  /* defining a mean "velocity" for static calculations: is used to
185  extrapolate the present results for next increment */
186 
187  if(*nmethod != 4){
188  for(i=0;i<*nk;i++){
189  for(j=1;j<mt;j++){
190  veold[mt*i+j]=(vold[mt*i+j]-vini[mt*i+j])/(*dtime);
191  }
192  }
193  }
194 
195  /* check whether next increment size must be decreased */
196 
197  if((*iit>il)&&(*idrct==0)){
198  if(*mortar==0){
199  *dtheta=*dthetaref*db;
200  *dthetaref=*dtheta;
201  printf(" convergence; the increment size is decreased to %e\n\n",*dtheta**tper);
202  if(*dtheta<*tmin){
203  printf("\n *ERROR: increment size smaller than minimum\n");
204  printf(" best solution and residuals are in the frd file\n\n");
205  NNEW(fn,double,mt**nk);
206  NNEW(inum,ITG,*nk);for(k=0;k<*nk;k++) inum[k]=1;
207  FORTRAN(storeresidual,(nactdof,b,fn,filab,ithermal,
208  nk,sti,stn,ipkon,inum,kon,lakon,ne,mi,orab,
209  ielorien,co,itg,ntg,vold,ielmat,thicke,ielprop,prop));
210  ++*kode;
211 
212  (*ttime)+=(*time);
213  frd(co,nk,kon,ipkon,lakon,ne,vold,stn,inum,nmethod,
214  kode,filab,een,t1act,fn,ttime,epn,ielmat,matname,enern,
215  xstaten,nstate_,istep,iinc,ithermal,qfn,mode,noddiam,
216  trab,inotr,ntrans,orab,ielorien,norien,description,
217  ipneigh,neigh,mi,sti,vr,vi,stnr,stni,vmax,stnmax,
218  &ngraph,veold,ener,ne,cs,set,nset,istartset,iendset,
219  ialset,eenmax,fnr,fni,emn,thicke,jobnamec,output,qfx,
220  cdn,mortar,cdnr,cdni,nmat);
221 
222  FORTRAN(stop,());
223  }
224  }
225  else{
226  printf("convergence\n\n");}
227  }
228 
229  /* check whether next increment size can be increased */
230 
231  else if(*iit<=ig){
232  if((*istab==1)&&(*idrct==0)){
233  *dtheta=*dthetaref*dd;
234  *dthetaref=*dtheta;
235  printf(" convergence; the increment size is increased to %e\n\n",*dtheta**tper);
236  }
237  else{
238  *istab=1;
239  printf(" convergence\n\n");
240  *dtheta=*dthetaref;
241  }
242  }
243  else{
244  *istab=0;
245  printf(" convergence\n\n");
246  *dtheta=*dthetaref;
247  }
248 
249  /* check whether new increment size exceeds maximum increment
250  size allowed (by the user) */
251 
252  if((*dtheta>*tmax)&&(*idrct==0)){
253  *dtheta=*tmax;
254  *dthetaref=*dtheta;
255  printf(" the increment size exceeds thetamax and is decreased to %e\n\n",*dtheta**tper);
256  }
257 
258  /* if itp=1 the increment just finished ends at a time point */
259 
260  if((*itpamp>0)&&(*idrct==0)){
261  if(*itp==1){
262  *jprint=*jout;
263  }else{
264  *jprint=*jout+1;
265  }
266  if(namta[3**itpamp-1]<0){
267  reftime=*ttime+*time+*dtheta**tper;
268  }else{
269  reftime=*time+*dtheta**tper;
270  }
271  istart=namta[3**itpamp-3];
272  iend=namta[3**itpamp-2];
273  FORTRAN(identamta,(amta,&reftime,&istart,&iend,&id));
274  if(id<istart){
275  inew=istart;
276  }else{
277  inew=id+1;
278  }
279 
280  /* inew: smallest time point exceeding time+dtheta*tper
281  inext: smallest time point exceeding time */
282 
283  if(*inext<inew){
284  if(namta[3**itpamp-1]<0){
285  *dtheta=(amta[2**inext-2]-*ttime-*time)/(*tper);
286  }else{
287  *dtheta=(amta[2**inext-2]-*time)/(*tper);
288  }
289  (*inext)++;
290  *itp=1;
291  printf(" the increment size exceeds a time point and is decreased to %e\n\n",*dtheta**tper);
292  }else{*itp=0;}
293  }
294 
295  /* check whether new time point exceeds end of step */
296 
297  if(*dtheta>=1.-*theta){
298  if(*dtheta>1.-*theta){iexceed=1;}else{iexceed=0;}
299 // if((*mortar==1)&&(1.-*theta>0.01)){
300 // *dtheta=0.999-*theta;
301 // if((*mortar==1)&&(1.-*theta>0.5)){
302 // *dtheta=0.900-*theta;
303 // }else{
304  *dtheta=1.-*theta;
305 // }
306  *dthetaref=*dtheta;
307  if(iexceed==1)
308  printf(" the increment size exceeds the remainder of the step and is decreased to %e\n\n",*dtheta**tper);
309 // if(*dtheta<=1.e-6){
310 // (*ttime)+=(*dtheta**tper);
311 // *dtime=0.;
312 // }
313  }
314  }
315  else{
316 
317  /* no convergence */
318 
319  /* check for the amount of iterations */
320 
321  if((*iit>ic)&&(*mortar==0)){
322  printf("\n *ERROR: too many iterations needed\n");
323  printf(" best solution and residuals are in the frd file\n\n");
324 
325  FORTRAN(writesummarydiv,(istep,iinc,icutb,iit,ttime,time,dtime));
326 
327  NNEW(fn,double,mt**nk);
328  NNEW(inum,ITG,*nk);for(k=0;k<*nk;k++) inum[k]=1;
329  FORTRAN(storeresidual,(nactdof,b,fn,filab,ithermal,nk,sti,stn,
330  ipkon,inum,kon,lakon,ne,mi,orab,ielorien,co,itg,ntg,vold,
331  ielmat,thicke,ielprop,prop));
332  ++*kode;
333 
334  (*ttime)+=(*time);
335  frd(co,nk,kon,ipkon,lakon,ne,vold,stn,inum,nmethod,
336  kode,filab,een,t1act,fn,ttime,epn,ielmat,matname,enern,
337  xstaten,nstate_,istep,iinc,ithermal,qfn,mode,noddiam,
338  trab,inotr,ntrans,orab,ielorien,norien,description,
339  ipneigh,neigh,mi,sti,vr,vi,stnr,stni,vmax,stnmax,
340  &ngraph,veold,ener,ne,cs,set,nset,istartset,iendset,
341  ialset,eenmax,fnr,fni,emn,thicke,jobnamec,output,qfx,cdn,
342  mortar,cdnr,cdni,nmat);
343 
344  FORTRAN(stop,());
345  }
346 
347  /* check for diverging residuals */
348 
349  if((*iit>=i0)||(fabs(ram[0])>1.e20)||(fabs(cam[0])>1.e20)||
350  (fabs(ram[1])>1.e20)||(fabs(cam[1])>1.e20)||
351  (cam[2]>*deltmx)||(qa[2]>0.)){
352  if((*ithermal!=2)&&(*mortar==0)){
353  if((ram1[0]>ram2[0])&&(ram[0]>ram2[0])&&(ram[0]>c1[0]*qam[0]))
354  idivergence=1;
355  }
356 
357 /* if((*ithermal!=2)&&(*mortar==1)){
358  if((ram[0]>ram[6])&&(ram1[0]>ram[6])){
359  printf("divergence allowed: residual force not converging\n");
360  if((ram1[0]>ram2[0])&&(ram[0]>ram2[0])&&(ram[0]>c1[0]*qam[0]))
361  idivergence=1;
362  }
363  if(ram[5]<0.5){
364  printf("divergence allowed: number of contact elements stabilized\n");
365  if((ram1[0]>ram2[0])&&(ram[0]>ram2[0])&&(ram[0]>c1[0]*qam[0]))
366  idivergence=1;
367  }
368  if(ram[5]>=ram1[3]){
369  if(((ram[4]>0.97*ram1[2])&&(ram[4]<1.03*ram1[2]))&&
370  ((ram[4]>0.97*ram2[2])&&(ram[4]<1.03*ram2[2]))){
371  printf("divergence allowed: repetitive pattern detected\n");
372  if((ram1[0]>ram2[0])&&(ram[0]>ram2[0])&&(ram[0]>c1[0]*qam[0]))
373  idivergence=1;
374  if(((ram[4]>0.99*ram1[2])&&(ram[4]<1.01*ram1[2]))&&
375  ((ram[4]>0.99*ram2[2])&&(ram[4]<1.01*ram2[2])&&(ram[0]>c1[0]*qam[0]))){
376  printf("repetitive pattern within 1 per cent -> divergence\n");
377  idivergence=1;
378  }
379  if((ram[4]>(1-1.e-3)*ram1[2])&&(ram[4]<(1+1.e-3)*ram1[2])&&
380  (ram[4]>(1-1.e-3)*ram2[2])&&(ram[4]<(1+1.e-3)*ram2[2])&&
381  (ram[5]==ram1[3])&&(ram[5]==ram2[3])){
382  printf("infinite loop -> divergence\n");
383  idivergence=1;
384  }
385  }
386  }
387  }*/
388 
389  if((*ithermal!=2)&&(*mortar==1)){
390 
391  if(ram[0]>1.e9){
392  printf("divergence allowed: residual force too large\n");
393  if((ram1[0]>ram2[0])&&(ram[0]>ram2[0])&&(ram[0]>c1[0]*qam[0]))
394  idivergence=1;
395  }
396 
397  /* number of contact elements does not change */
398 
399  if(((ITG)ram[6]==0)&&((ITG)ram1[6]==0)){
400  printf("divergence allowed: number of contact elements stabilized\n");
401  if((ram1[0]>ram2[0])&&(ram[0]>ram2[0])&&(ram[0]>c1[0]*qam[0]))
402  idivergence=1;
403  }
404 
405  /* rate of number of contact elements is increasing */
406 
407  if(((ITG)ram[6]>=(ITG)ram1[6])&&((ITG)ram[6]>=(ITG)ram2[6])){
408 
409  if(((ram[4]>0.98*ram1[4])&&(ram[4]<1.02*ram1[4]))&&
410  ((ram[4]>0.98*ram2[4])&&(ram[4]<1.02*ram2[4]))){
411  printf("divergence allowed: repetitive pattern detected\n");
412  if((ram1[0]>ram2[0])&&(ram[0]>ram2[0])&&(ram[0]>c1[0]*qam[0]))
413  idivergence=1;
414  }
415 
416  /* if(((ram[4]>0.99*ram1[4])&&(ram[4]<1.01*ram1[4]))&&
417  ((ram[4]>0.99*ram2[4])&&(ram[4]<1.01*ram2[4]))&&
418  (ram[0]>c1[0]*qam[0])){
419  printf("repetitive pattern within 1 %: divergence\n");
420  idivergence=1;
421  }
422 
423  if(((ram[4]>(1.-1.e-3)*ram1[4])&&(ram[4]<(1.+1.e-3)*ram1[4]))&&
424  ((ram[4]>(1.-1.e-3)*ram2[4])&&(ram[4]<(1.+1.e-3)*ram2[4]))&&
425  ((ITG)ram[6]==(ITG)ram1[6])&&((ITG)ram[6]==(ITG)ram2[6])){
426  printf("infinite loop -> divergence\n");
427  idivergence=1;
428  }*/
429  }
430 
431  /* if(((ITG)ram1[7]==1)&&((ITG)ram2[7]!=-1)){
432  printf("divergence if ram>ram1>ram2: contact elements not stabilizing\n");
433  if((ram1[0]>ram2[0])&&(ram[0]>ram2[0])&&(ram[0]>c1[0]*qam[0]))
434  idivergence=1;
435  }*/
436  }
437 
438 
439  /* for thermal calculations the maximum temperature change
440  is checked as well */
441 
442  if(*ithermal>1){
443  if((ram1[1]>ram2[1])&&(ram[1]>ram2[1])&&(ram[1]>c1[1]*qam[1]))
444  idivergence=1;
445  if(cam[2]>*deltmx) idivergence=2;
446  }
447  if(idivergence>0){
448  if(*idrct==1) {
449 
450  /* fixed time increments */
451 
452  printf("\n *ERROR: solution seems to diverge; please try \n");
453  printf(" automatic incrementation; program stops\n");
454  printf(" best solution and residuals are in the frd file\n\n");
455 
456  FORTRAN(writesummarydiv,(istep,iinc,icutb,iit,ttime,time,
457  dtime));
458 
459  NNEW(fn,double,mt**nk);
460  NNEW(inum,ITG,*nk);for(k=0;k<*nk;k++) inum[k]=1;
461  FORTRAN(storeresidual,(nactdof,b,fn,filab,ithermal,nk,
462  sti,stn,ipkon,inum,kon,lakon,ne,mi,orab,
463  ielorien,co,itg,ntg,vold,ielmat,thicke,ielprop,prop));
464  ++*kode;
465 
466  (*ttime)+=(*time);
467  frd(co,nk,kon,ipkon,lakon,ne,vold,stn,inum,nmethod,
468  kode,filab,een,t1act,fn,ttime,epn,ielmat,matname,enern,
469  xstaten,nstate_,istep,iinc,ithermal,qfn,mode,noddiam,
470  trab,inotr,ntrans,orab,ielorien,norien,description,
471  ipneigh,neigh,mi,sti,vr,vi,stnr,stni,vmax,stnmax,
472  &ngraph,veold,ener,ne,cs,set,nset,istartset,iendset,
473  ialset,eenmax,fnr,fni,emn,thicke,jobnamec,output,qfx,cdn,
474  mortar,cdnr,cdni,nmat);
475 
476  FORTRAN(stop,());
477  }
478  else {
479 
480  /* variable time increments */
481 
482  if(qa[2]>0.){
483  *dtheta=*dtheta*qa[2];
484  printf("increment size decrease requested by a material user routine (through pnewdt)\n\n");
485  }else{
486  if(idivergence==1){
487  *dtheta=*dtheta*df;
488  }else{
489  *dtheta=*dtheta**deltmx/cam[2]*da;
490  }
491  }
492  *dthetaref=*dtheta;
493  printf(" divergence; the increment size is decreased to %e\n",*dtheta**tper);
494  printf(" the increment is reattempted\n\n");
495 
496  FORTRAN(writesummarydiv,(istep,iinc,icutb,iit,ttime,time,
497  dtime));
498 
499  *istab=0;
500  if(*itp==1){
501  *itp=0;
502  (*inext)--;
503  }
504 
505  /* check whether new increment size is smaller than minimum */
506 
507  if(*dtheta<*tmin){
508  printf("\n *ERROR: increment size smaller than minimum\n");
509  printf(" best solution and residuals are in the frd file\n\n");
510  NNEW(fn,double,mt**nk);
511  NNEW(inum,ITG,*nk);for(k=0;k<*nk;k++) inum[k]=1;
512  FORTRAN(storeresidual,(nactdof,b,fn,filab,ithermal,
513  nk,sti,stn,ipkon,inum,kon,lakon,ne,mi,orab,
514  ielorien,co,itg,ntg,vold,ielmat,thicke,ielprop,prop));
515  ++*kode;
516 
517  (*ttime)+=(*time);
518  frd(co,nk,kon,ipkon,lakon,ne,vold,stn,inum,nmethod,
519  kode,filab,een,t1act,fn,ttime,epn,ielmat,matname,enern,
520  xstaten,nstate_,istep,iinc,ithermal,qfn,mode,noddiam,
521  trab,inotr,ntrans,orab,ielorien,norien,description,
522  ipneigh,neigh,mi,sti,vr,vi,stnr,stni,vmax,stnmax,
523  &ngraph,veold,ener,ne,cs,set,nset,istartset,iendset,
524  ialset,eenmax,fnr,fni,emn,thicke,jobnamec,output,qfx,cdn,
525  mortar,cdnr,cdni,nmat);
526 
527  FORTRAN(stop,());
528  }
529  *icntrl=1;
530  (*icutb)++;
531 
532  /* check whether too many cutbacks */
533 
534  if(*icutb>ia){
535  printf("\n *ERROR: too many cutbacks\n");
536  printf(" best solution and residuals are in the frd file\n\n");
537  NNEW(fn,double,mt**nk);
538  NNEW(inum,ITG,*nk);for(k=0;k<*nk;k++) inum[k]=1;
539  FORTRAN(storeresidual,(nactdof,b,fn,filab,ithermal,
540  nk,sti,stn,ipkon,inum,kon,lakon,ne,mi,orab,
541  ielorien,co,itg,ntg,vold,ielmat,thicke,ielprop,prop));
542  ++*kode;
543 
544  (*ttime)+=(*time);
545  frd(co,nk,kon,ipkon,lakon,ne,vold,stn,inum,nmethod,
546  kode,filab,een,t1act,fn,ttime,epn,ielmat,matname,enern,
547  xstaten,nstate_,istep,iinc,ithermal,qfn,mode,noddiam,
548  trab,inotr,ntrans,orab,ielorien,norien,description,
549  ipneigh,neigh,mi,sti,vr,vi,stnr,stni,vmax,stnmax,
550  &ngraph,veold,ener,ne,cs,set,nset,istartset,iendset,
551  ialset,eenmax,fnr,fni,emn,thicke,jobnamec,output,qfx,cdn,
552  mortar,cdnr,cdni,nmat);
553 
554  FORTRAN(stop,());
555  }
556  if(*uncoupled){
557  if(*ithermal==1){
558  (ctrl[0])/=4;
559  }
560  *ithermal=3;
561  }
562 
563  /* default value for qa[2] */
564 
565  qa[2]=-1.;
566 
567  return;
568  }
569  }
570  }
571 
572  /* check for too slow convergence */
573 
574  if(*iit>=ir){
575  if(*ithermal!=2){
576  iest1=(ITG)ceil(*iit+log(ran*qam[0]/(ram[0]))/
577  log(ram[0]/(ram1[0])));
578  }
579  if(*ithermal>1){
580  iest2=(ITG)ceil(*iit+log(ran*qam[1]/(ram[1]))/
581  log(ram[1]/(ram1[1])));
582  }
583  if(iest1>iest2){iest=iest1;}else{iest=iest2;}
584  if((iest>0)&&(*mortar==0)){
585  printf(" estimated number of iterations till convergence = %" ITGFORMAT "\n",
586  iest);
587  }
588  if((((iest>ic)||(*iit==ic))&&(*mortar==0))||((*mortar==1)&&(*iit==60))){
589 
590  if(*idrct!=1){
591  *dtheta=*dtheta*dc;
592  *dthetaref=*dtheta;
593  printf(" too slow convergence; the increment size is decreased to %e\n",*dtheta**tper);
594  printf(" the increment is reattempted\n\n");
595 
596  FORTRAN(writesummarydiv,(istep,iinc,icutb,iit,ttime,
597  time,dtime));
598 
599  *istab=0;
600 
601  /* check whether new increment size is smaller than minimum */
602 
603  if(*dtheta<*tmin){
604  printf("\n *ERROR: increment size smaller than minimum\n");
605  printf(" best solution and residuals are in the frd file\n\n");
606  NNEW(fn,double,mt**nk);
607  NNEW(inum,ITG,*nk);for(k=0;k<*nk;k++) inum[k]=1;
608  FORTRAN(storeresidual,(nactdof,b,fn,filab,ithermal,
609  nk,sti,stn,ipkon,inum,kon,lakon,ne,mi,orab,
610  ielorien,co,itg,ntg,vold,ielmat,thicke,ielprop,prop));
611  ++*kode;
612 
613  (*ttime)+=(*time);
614  frd(co,nk,kon,ipkon,lakon,ne,vold,stn,inum,nmethod,
615  kode,filab,een,t1act,fn,ttime,epn,ielmat,matname,enern,
616  xstaten,nstate_,istep,iinc,ithermal,qfn,mode,noddiam,
617  trab,inotr,ntrans,orab,ielorien,norien,description,
618  ipneigh,neigh,mi,sti,vr,vi,stnr,stni,vmax,stnmax,
619  &ngraph,veold,ener,ne,cs,set,nset,istartset,iendset,
620  ialset,eenmax,fnr,fni,emn,thicke,jobnamec,output,qfx,cdn,
621  mortar,cdnr,cdni,nmat);
622 
623  FORTRAN(stop,());
624  }
625  *icntrl=1;
626  (*icutb)++;
627  if(*icutb>ia){
628  printf("\n *ERROR: too many cutbacks\n");
629  printf(" best solution and residuals are in the frd file\n\n");
630  NNEW(fn,double,mt**nk);
631  NNEW(inum,ITG,*nk);for(k=0;k<*nk;k++) inum[k]=1;
632  FORTRAN(storeresidual,(nactdof,b,fn,filab,ithermal,
633  nk,sti,stn,ipkon,inum,kon,lakon,ne,mi,orab,
634  ielorien,co,itg,ntg,vold,ielmat,thicke,ielprop,prop));
635  ++*kode;
636 
637  (*ttime)+=(*time);
638  frd(co,nk,kon,ipkon,lakon,ne,vold,stn,inum,nmethod,
639  kode,filab,een,t1act,fn,ttime,epn,ielmat,matname,enern,
640  xstaten,nstate_,istep,iinc,ithermal,qfn,mode,noddiam,
641  trab,inotr,ntrans,orab,ielorien,norien,description,
642  ipneigh,neigh,mi,sti,vr,vi,stnr,stni,vmax,stnmax,
643  &ngraph,veold,ener,ne,cs,set,nset,istartset,iendset,
644  ialset,eenmax,fnr,fni,emn,thicke,jobnamec,output,qfx,cdn,
645  mortar,cdnr,cdni,nmat);
646 
647  FORTRAN(stop,());
648  }
649  if(*uncoupled){
650  if(*ithermal==1){
651  (ctrl[0])/=4;
652  }
653  *ithermal=3;
654  }
655  return;
656  }
657  }
658  }
659 
660  printf(" no convergence\n\n");
661 
662  (*iit)++;
663 
664  }
665 
666  /* default value for qa[2] */
667 
668  /* qa[2]=-1;*/
669 
670  return;
671 }
#define ITGFORMAT
Definition: CalculiX.h:52
void frd(double *co, ITG *nk, ITG *kon, ITG *ipkon, char *lakon, ITG *ne0, double *v, double *stn, ITG *inum, ITG *nmethod, ITG *kode, char *filab, double *een, double *t1, double *fn, double *time, double *epn, ITG *ielmat, char *matname, double *enern, double *xstaten, ITG *nstate_, ITG *istep, ITG *iinc, ITG *ithermal, double *qfn, ITG *mode, ITG *noddiam, double *trab, ITG *inotr, ITG *ntrans, double *orab, ITG *ielorien, ITG *norien, char *description, ITG *ipneigh, ITG *neigh, ITG *mi, double *stx, double *vr, double *vi, double *stnr, double *stni, double *vmax, double *stnmax, ITG *ngraph, double *veold, double *ener, ITG *ne, double *cs, char *set, ITG *nset, ITG *istartset, ITG *iendset, ITG *ialset, double *eenmax, double *fnr, double *fni, double *emn, double *thicke, char *jobnamec, char *output, double *qfx, double *cdn, ITG *mortar, double *cdnr, double *cdni, ITG *nmat)
Definition: frd.c:27
void FORTRAN(addimdnodecload,(ITG *nodeforc, ITG *i, ITG *imdnode, ITG *nmdnode, double *xforc, ITG *ikmpc, ITG *ilmpc, ITG *ipompc, ITG *nodempc, ITG *nmpc, ITG *imddof, ITG *nmddof, ITG *nactdof, ITG *mi, ITG *imdmpc, ITG *nmdmpc, ITG *imdboun, ITG *nmdboun, ITG *ikboun, ITG *nboun, ITG *ilboun, ITG *ithermal))
subroutine df(x, u, uprime, rpar, nev)
Definition: subspace.f:132
subroutine writesummary(istep, j, icutb, l, ttime, time, dtime)
Definition: writesummary.f:19
subroutine storeresidual(nactdof, b, fn, filab, ithermal, nk, sti, stn, ipkon, inum, kon, lakon, ne, mi, orab, ielorien, co, itg, ntg, vold, ielmat, thicke, ielprop, prop)
Definition: storeresidual.f:19
subroutine stop()
Definition: stop.f:19
subroutine writesummarydiv(istep, j, icutb, l, ttime, time, dtime)
Definition: writesummarydiv.f:19
subroutine identamta(amta, reftime, istart, iend, id)
Definition: identamta.f:25
#define ITG
Definition: CalculiX.h:51
#define NNEW(a, b, c)
Definition: CalculiX.h:39
Hosted by OpenAircraft.com, (Michigan UAV, LLC)