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

Go to the source code of this file.

Functions

void linstatic (double *co, ITG *nk, ITG *kon, ITG *ipkon, char *lakon, ITG *ne, ITG *nodeboun, ITG *ndirboun, double *xboun, ITG *nboun, ITG *ipompc, ITG *nodempc, double *coefmpc, char *labmpc, ITG *nmpc, ITG *nodeforc, ITG *ndirforc, double *xforc, ITG *nforc, ITG *nelemload, char *sideload, double *xload, ITG *nload, ITG *nactdof, ITG **icolp, ITG *jq, ITG **irowp, ITG *neq, ITG *nzl, ITG *nmethod, ITG *ikmpc, ITG *ilmpc, ITG *ikboun, ITG *ilboun, double *elcon, ITG *nelcon, double *rhcon, ITG *nrhcon, double *alcon, ITG *nalcon, double *alzero, ITG *ielmat, ITG *ielorien, ITG *norien, double *orab, ITG *ntmat_, double *t0, double *t1, double *t1old, ITG *ithermal, double *prestr, ITG *iprestr, double *vold, ITG *iperturb, double *sti, ITG *nzs, ITG *kode, char *filab, double *eme, ITG *iexpl, double *plicon, ITG *nplicon, double *plkcon, ITG *nplkcon, double *xstate, ITG *npmat_, char *matname, ITG *isolver, ITG *mi, ITG *ncmat_, ITG *nstate_, double *cs, ITG *mcs, ITG *nkon, double *ener, double *xbounold, double *xforcold, double *xloadold, char *amname, double *amta, ITG *namta, ITG *nam, ITG *iamforc, ITG *iamload, ITG *iamt1, ITG *iamboun, double *ttime, char *output, char *set, ITG *nset, ITG *istartset, ITG *iendset, ITG *ialset, ITG *nprint, char *prlab, char *prset, ITG *nener, double *trab, ITG *inotr, ITG *ntrans, double *fmpc, char *cbody, ITG *ibody, double *xbody, ITG *nbody, double *xbodyold, double *timepar, double *thicke, char *jobnamec, char *tieset, ITG *ntie, ITG *istep, ITG *nmat, ITG *ielprop, double *prop, char *typeboun)
 

Function Documentation

void linstatic ( double *  co,
ITG nk,
ITG kon,
ITG ipkon,
char *  lakon,
ITG ne,
ITG nodeboun,
ITG ndirboun,
double *  xboun,
ITG nboun,
ITG ipompc,
ITG nodempc,
double *  coefmpc,
char *  labmpc,
ITG nmpc,
ITG nodeforc,
ITG ndirforc,
double *  xforc,
ITG nforc,
ITG nelemload,
char *  sideload,
double *  xload,
ITG nload,
ITG nactdof,
ITG **  icolp,
ITG jq,
ITG **  irowp,
ITG neq,
ITG nzl,
ITG nmethod,
ITG ikmpc,
ITG ilmpc,
ITG ikboun,
ITG ilboun,
double *  elcon,
ITG nelcon,
double *  rhcon,
ITG nrhcon,
double *  alcon,
ITG nalcon,
double *  alzero,
ITG ielmat,
ITG ielorien,
ITG norien,
double *  orab,
ITG ntmat_,
double *  t0,
double *  t1,
double *  t1old,
ITG ithermal,
double *  prestr,
ITG iprestr,
double *  vold,
ITG iperturb,
double *  sti,
ITG nzs,
ITG kode,
char *  filab,
double *  eme,
ITG iexpl,
double *  plicon,
ITG nplicon,
double *  plkcon,
ITG nplkcon,
double *  xstate,
ITG npmat_,
char *  matname,
ITG isolver,
ITG mi,
ITG ncmat_,
ITG nstate_,
double *  cs,
ITG mcs,
ITG nkon,
double *  ener,
double *  xbounold,
double *  xforcold,
double *  xloadold,
char *  amname,
double *  amta,
ITG namta,
ITG nam,
ITG iamforc,
ITG iamload,
ITG iamt1,
ITG iamboun,
double *  ttime,
char *  output,
char *  set,
ITG nset,
ITG istartset,
ITG iendset,
ITG ialset,
ITG nprint,
char *  prlab,
char *  prset,
ITG nener,
double *  trab,
ITG inotr,
ITG ntrans,
double *  fmpc,
char *  cbody,
ITG ibody,
double *  xbody,
ITG nbody,
double *  xbodyold,
double *  timepar,
double *  thicke,
char *  jobnamec,
char *  tieset,
ITG ntie,
ITG istep,
ITG nmat,
ITG ielprop,
double *  prop,
char *  typeboun 
)

Definition at line 35 of file linstatic.c.

68  {
69 
70  char description[13]=" ";
71 
72  ITG *inum=NULL,k,*icol=NULL,*irow=NULL,ielas,icmd=0,iinc=1,nasym=0,i,j,ic,ir,
73  mass[2]={0,0}, stiffness=1, buckling=0, rhsi=1, intscheme=0,*ncocon=NULL,
74  *nshcon=NULL,mode=-1,noddiam=-1,*ipobody=NULL,inewton=0,coriolis=0,iout,
75  ifreebody,*itg=NULL,ntg=0,symmetryflag=0,inputformat=0,ngraph=1,im,
76  mt=mi[1]+1,ne0,*integerglob=NULL,iglob=0,*ipneigh=NULL,*neigh=NULL,
77  icfd=0,*inomat=NULL,mortar,*islavact=NULL,*islavnode=NULL,*nslavnode=NULL,
78  *islavsurf=NULL,nretain,*iretain=NULL,*noderetain=NULL,*ndirretain=NULL,
79  nmethodl;
80 
81  double *stn=NULL,*v=NULL,*een=NULL,cam[5],*xstiff=NULL,*stiini=NULL,*tper,
82  *f=NULL,*fn=NULL,qa[3],*fext=NULL,*epn=NULL,*xstateini=NULL,
83  *vini=NULL,*stx=NULL,*enern=NULL,*xbounact=NULL,*xforcact=NULL,
84  *xloadact=NULL,*t1act=NULL,*ampli=NULL,*xstaten=NULL,*eei=NULL,
85  *enerini=NULL,*cocon=NULL,*shcon=NULL,*physcon=NULL,*qfx=NULL,
86  *qfn=NULL,sigma=0.,*cgr=NULL,*xbodyact=NULL,*vr=NULL,*vi=NULL,
87  *stnr=NULL,*stni=NULL,*vmax=NULL,*stnmax=NULL,*springarea=NULL,
88  *eenmax=NULL,*fnr=NULL,*fni=NULL,*emn=NULL,*clearini=NULL,ptime,
89  *emeini=NULL,*doubleglob=NULL,*au=NULL,*ad=NULL,*b=NULL,*aub=NULL,
90  *adb=NULL,*pslavsurf=NULL,*pmastsurf=NULL,*cdn=NULL,*cdnr=NULL,
91  *cdni=NULL,*submatrix=NULL;
92 
93 #ifdef SGI
94  ITG token;
95 #endif
96 
97  /* dummy arguments for the results call */
98 
99  double *veold=NULL,*accold=NULL,bet,gam,dtime,time,reltime=1.;
100 
101  icol=*icolp;
102  irow=*irowp;
103 
104  tper=&timepar[1];
105 
106  time=*tper;
107  dtime=*tper;
108 
109  ne0=*ne;
110 
111  /* determining the global values to be used as boundary conditions
112  for a submodel */
113 
114  getglobalresults(jobnamec,&integerglob,&doubleglob,nboun,iamboun,xboun,
115  nload,sideload,iamload,&iglob,nforc,iamforc,xforc,
116  ithermal,nk,t1,iamt1);
117 
118  /* allocating fields for the actual external loading */
119 
120  NNEW(xbounact,double,*nboun);
121  for(k=0;k<*nboun;++k){xbounact[k]=xbounold[k];}
122  NNEW(xforcact,double,*nforc);
123  NNEW(xloadact,double,2**nload);
124  NNEW(xbodyact,double,7**nbody);
125  /* copying the rotation axis and/or acceleration vector */
126  for(k=0;k<7**nbody;k++){xbodyact[k]=xbody[k];}
127  if(*ithermal==1){
128  NNEW(t1act,double,*nk);
129  for(k=0;k<*nk;++k){t1act[k]=t1old[k];}
130  }
131 
132  /* assigning the body forces to the elements */
133 
134  if(*nbody>0){
135  ifreebody=*ne+1;
136  NNEW(ipobody,ITG,2*ifreebody**nbody);
137  for(k=1;k<=*nbody;k++){
138  FORTRAN(bodyforce,(cbody,ibody,ipobody,nbody,set,istartset,
139  iendset,ialset,&inewton,nset,&ifreebody,&k));
140  RENEW(ipobody,ITG,2*(*ne+ifreebody));
141  }
142  RENEW(ipobody,ITG,2*(ifreebody-1));
143  }
144 
145  /* allocating a field for the instantaneous amplitude */
146 
147  NNEW(ampli,double,*nam);
148 
149  FORTRAN(tempload,(xforcold,xforc,xforcact,iamforc,nforc,xloadold,xload,
150  xloadact,iamload,nload,ibody,xbody,nbody,xbodyold,xbodyact,
151  t1old,t1,t1act,iamt1,nk,amta,
152  namta,nam,ampli,&time,&reltime,ttime,&dtime,ithermal,nmethod,
153  xbounold,xboun,xbounact,iamboun,nboun,
154  nodeboun,ndirboun,nodeforc,ndirforc,istep,&iinc,
155  co,vold,itg,&ntg,amname,ikboun,ilboun,nelemload,sideload,mi,
156  ntrans,trab,inotr,veold,integerglob,doubleglob,tieset,istartset,
157  iendset,ialset,ntie,nmpc,ipompc,ikmpc,ilmpc,nodempc,coefmpc));
158 
159  /* determining the internal forces and the stiffness coefficients */
160 
161  NNEW(f,double,*neq);
162 
163  /* allocating a field for the stiffness matrix */
164 
165  NNEW(xstiff,double,(long long)27*mi[0]**ne);
166 
167  iout=-1;
168  NNEW(v,double,mt**nk);
169  NNEW(fn,double,mt**nk);
170  NNEW(stx,double,6*mi[0]**ne);
171  NNEW(inum,ITG,*nk);
172  results(co,nk,kon,ipkon,lakon,ne,v,stn,inum,stx,
173  elcon,nelcon,rhcon,nrhcon,alcon,nalcon,alzero,ielmat,
174  ielorien,norien,orab,ntmat_,t0,t1act,ithermal,
175  prestr,iprestr,filab,eme,emn,een,iperturb,
176  f,fn,nactdof,&iout,qa,vold,b,nodeboun,
177  ndirboun,xbounact,nboun,ipompc,
178  nodempc,coefmpc,labmpc,nmpc,nmethod,cam,neq,veold,accold,
179  &bet,&gam,&dtime,&time,ttime,plicon,nplicon,plkcon,nplkcon,
180  xstateini,xstiff,xstate,npmat_,epn,matname,mi,&ielas,
181  &icmd,ncmat_,nstate_,stiini,vini,ikboun,ilboun,ener,enern,
182  emeini,xstaten,eei,enerini,cocon,ncocon,set,nset,istartset,
183  iendset,ialset,nprint,prlab,prset,qfx,qfn,trab,inotr,ntrans,
184  fmpc,nelemload,nload,ikmpc,ilmpc,istep,&iinc,springarea,
185  &reltime,&ne0,xforc,nforc,thicke,shcon,nshcon,
186  sideload,xloadact,xloadold,&icfd,inomat,pslavsurf,pmastsurf,
187  &mortar,islavact,cdn,islavnode,nslavnode,ntie,clearini,
188  islavsurf,ielprop,prop);
189  SFREE(v);SFREE(fn);SFREE(stx);SFREE(inum);
190  iout=1;
191 
192  /* determining the system matrix and the external forces */
193 
194  NNEW(ad,double,*neq);
195  NNEW(fext,double,*neq);
196 
197  if(*nmethod==11){
198  NNEW(au,double,nzs[2]);
199  rhsi=0;
200  nmethodl=2;
201  }else{
202  NNEW(au,double,*nzs);
203  nmethodl=*nmethod;
204  }
205 
206 
207  FORTRAN(mafillsm,(co,nk,kon,ipkon,lakon,ne,nodeboun,ndirboun,xbounact,nboun,
208  ipompc,nodempc,coefmpc,nmpc,nodeforc,ndirforc,xforcact,
209  nforc,nelemload,sideload,xloadact,nload,xbodyact,ipobody,
210  nbody,cgr,ad,au,fext,nactdof,icol,jq,irow,neq,nzl,&nmethodl,
211  ikmpc,ilmpc,ikboun,ilboun,
212  elcon,nelcon,rhcon,nrhcon,alcon,nalcon,alzero,ielmat,
213  ielorien,norien,orab,ntmat_,
214  t0,t1act,ithermal,prestr,iprestr,vold,iperturb,sti,
215  nzs,stx,adb,aub,iexpl,plicon,nplicon,plkcon,nplkcon,
216  xstiff,npmat_,&dtime,matname,mi,
217  ncmat_,mass,&stiffness,&buckling,&rhsi,&intscheme,physcon,
218  shcon,nshcon,cocon,ncocon,ttime,&time,istep,&iinc,&coriolis,
219  ibody,xloadold,&reltime,veold,springarea,nstate_,
220  xstateini,xstate,thicke,integerglob,doubleglob,
221  tieset,istartset,iendset,ialset,ntie,&nasym,pslavsurf,
222  pmastsurf,&mortar,clearini,ielprop,prop));
223 
224  /* determining the right hand side */
225 
226  NNEW(b,double,*neq);
227  for(k=0;k<*neq;++k){
228  b[k]=fext[k]-f[k];
229  }
230  SFREE(fext);SFREE(f);
231 
232  /* generation of a substructure stiffness matrix */
233 
234  if(*nmethod==11){
235 
236  /* factorizing the matrix */
237 
238  if(*neq>0){
239  if(*isolver==0){
240 #ifdef SPOOLES
241  spooles_factor(ad,au,adb,aub,&sigma,icol,irow,neq,nzs,&symmetryflag,
242  &inputformat,&nzs[2]);
243 #else
244  printf("*ERROR in linstatic: the SPOOLES library is not linked\n\n");
245  FORTRAN(stop,());
246 #endif
247  }
248  else if(*isolver==7){
249 #ifdef PARDISO
250  pardiso_factor(ad,au,adb,aub,&sigma,icol,irow,neq,nzs,
251  &symmetryflag,&inputformat,jq,&nzs[2]);
252 #else
253  printf("*ERROR in linstatic: the PARDISO library is not linked\n\n");
254  FORTRAN(stop,());
255 #endif
256  }
257  }
258 
259  /* determining the nodes and the degrees of freedom in those nodes
260  belonging to the substructure */
261 
262  NNEW(iretain,ITG,*nk);
263  NNEW(noderetain,ITG,*nk);
264  NNEW(ndirretain,ITG,*nk);
265  nretain=0;
266 
267  for(i=0;i<*nboun;i++){
268  if(strcmp1(&typeboun[i],"C")==0){
269  iretain[nretain]=i+1;
270  noderetain[nretain]=nodeboun[i];
271  ndirretain[nretain]=ndirboun[i];
272  nretain++;
273  }
274  }
275 
276  RENEW(iretain,ITG,nretain);
277  RENEW(noderetain,ITG,nretain);
278  RENEW(ndirretain,ITG,nretain);
279 
280  /* solving the system of equations with appropriate rhs */
281 
282  NNEW(submatrix,double,nretain*nretain);
283 
284  for(i=0;i<nretain;i++){
285  DMEMSET(b,0,*neq,0.);
286  ic=*neq+iretain[i]-1;
287  for(j=jq[ic]-1;j<jq[ic+1]-1;j++){
288  ir=irow[j]-1;
289  b[ir]-=au[j];
290  }
291 
292  /* solving the system */
293 
294  if(*neq>0){
295  if(*isolver==0){
296 #ifdef SPOOLES
297  spooles_solve(b,neq);
298 #endif
299  }
300  else if(*isolver==7){
301 #ifdef PARDISO
302  pardiso_solve(b,neq,&symmetryflag);
303 #endif
304 
305  }
306  }
307 
308  /* calculating the internal forces */
309 
310  NNEW(v,double,mt**nk);
311  NNEW(fn,double,mt**nk);
312  NNEW(stn,double,6**nk);
313  NNEW(inum,ITG,*nk);
314  NNEW(stx,double,6*mi[0]**ne);
315 
316  if(strcmp1(&filab[261],"E ")==0) NNEW(een,double,6**nk);
317  if(strcmp1(&filab[522],"ENER")==0) NNEW(enern,double,*nk);
318  if(strcmp1(&filab[2697],"ME ")==0) NNEW(emn,double,6**nk);
319 
320  NNEW(eei,double,6*mi[0]**ne);
321  if(*nener==1){
322  NNEW(stiini,double,6*mi[0]**ne);
323  NNEW(enerini,double,mi[0]**ne);}
324 
325  /* replacing the appropriate boundary value by unity */
326 
327  xbounact[iretain[i]-1]=1.;
328 
329  results(co,nk,kon,ipkon,lakon,ne,v,stn,inum,stx,
330  elcon,nelcon,rhcon,nrhcon,alcon,nalcon,alzero,ielmat,
331  ielorien,norien,orab,ntmat_,t0,t1act,ithermal,
332  prestr,iprestr,filab,eme,emn,een,iperturb,
333  f,fn,nactdof,&iout,qa,vold,b,nodeboun,ndirboun,xbounact,nboun,ipompc,
334  nodempc,coefmpc,labmpc,nmpc,nmethod,cam,neq,veold,accold,&bet,
335  &gam,&dtime,&time,ttime,plicon,nplicon,plkcon,nplkcon,
336  xstateini,xstiff,xstate,npmat_,epn,matname,mi,&ielas,&icmd,
337  ncmat_,nstate_,stiini,vini,ikboun,ilboun,ener,enern,emeini,
338  xstaten,eei,enerini,cocon,ncocon,set,nset,istartset,iendset,
339  ialset,nprint,prlab,prset,qfx,qfn,trab,inotr,ntrans,fmpc,
340  nelemload,nload,ikmpc,ilmpc,istep,&iinc,springarea,&reltime,
341  &ne0,xforc,nforc,thicke,shcon,nshcon,
342  sideload,xloadact,xloadold,&icfd,inomat,pslavsurf,pmastsurf,
343  &mortar,islavact,cdn,islavnode,nslavnode,ntie,clearini,
344  islavsurf,ielprop,prop);
345 
346  xbounact[iretain[i]-1]=0.;
347 
348  SFREE(v);SFREE(stn);SFREE(inum);SFREE(stx);
349 
350  if(strcmp1(&filab[261],"E ")==0) SFREE(een);
351  if(strcmp1(&filab[522],"ENER")==0) SFREE(enern);
352  if(strcmp1(&filab[2697],"ME ")==0) SFREE(emn);
353 
354  SFREE(eei);if(*nener==1){SFREE(stiini);SFREE(enerini);}
355 
356  /* storing the internal forces in the substructure
357  stiffness matrix */
358 
359  for(j=0;j<nretain;j++){
360  submatrix[i*nretain+j]=fn[mt*(noderetain[j]-1)+ndirretain[j]];
361  }
362 
363  SFREE(fn);
364 
365  }
366 
367  SFREE(au);SFREE(ad);SFREE(b);
368  SFREE(iretain);
369 
370  SFREE(xbounact);SFREE(xforcact);SFREE(xloadact);SFREE(t1act);SFREE(ampli);
371  SFREE(xbodyact);if(*nbody>0) SFREE(ipobody);SFREE(xstiff);
372 
373  if(iglob==1){SFREE(integerglob);SFREE(doubleglob);}
374 
375  FORTRAN(writesubmatrix,(submatrix,noderetain,ndirretain,&nretain,jobnamec));
376 
377  SFREE(submatrix);SFREE(noderetain);SFREE(ndirretain);
378 
379  return;
380 
381 
382  }else if(*nmethod!=0){
383 
384  if(*isolver==0){
385 #ifdef SPOOLES
386  spooles(ad,au,adb,aub,&sigma,b,icol,irow,neq,nzs,&symmetryflag,
387  &inputformat,&nzs[2]);
388 #else
389  printf("*ERROR in linstatic: the SPOOLES library is not linked\n\n");
390  FORTRAN(stop,());
391 #endif
392  }
393  else if((*isolver==2)||(*isolver==3)){
394  preiter(ad,&au,b,&icol,&irow,neq,nzs,isolver,iperturb);
395  }
396  else if(*isolver==4){
397 #ifdef SGI
398  token=1;
399  sgi_main(ad,au,adb,aub,&sigma,b,icol,irow,neq,nzs,token);
400 #else
401  printf("*ERROR in linstatic: the SGI library is not linked\n\n");
402  FORTRAN(stop,());
403 #endif
404  }
405  else if(*isolver==5){
406 #ifdef TAUCS
407  tau(ad,&au,adb,aub,&sigma,b,icol,&irow,neq,nzs);
408 #else
409  printf("*ERROR in linstatic: the TAUCS library is not linked\n\n");
410  FORTRAN(stop,());
411 #endif
412  }
413  else if(*isolver==7){
414 #ifdef PARDISO
415  pardiso_main(ad,au,adb,aub,&sigma,b,icol,irow,neq,nzs,
416  &symmetryflag,&inputformat,jq,&nzs[2]);
417 #else
418  printf("*ERROR in linstatic: the PARDISO library is not linked\n\n");
419  FORTRAN(stop,());
420 #endif
421  }
422 
423  SFREE(ad);SFREE(au);
424 
425  /* calculating the displacements and the stresses and storing */
426  /* the results in frd format for each valid eigenmode */
427 
428  NNEW(v,double,mt**nk);
429  NNEW(fn,double,mt**nk);
430  NNEW(stn,double,6**nk);
431  NNEW(inum,ITG,*nk);
432  NNEW(stx,double,6*mi[0]**ne);
433 
434  if(strcmp1(&filab[261],"E ")==0) NNEW(een,double,6**nk);
435  if(strcmp1(&filab[522],"ENER")==0) NNEW(enern,double,*nk);
436  if(strcmp1(&filab[2697],"ME ")==0) NNEW(emn,double,6**nk);
437 
438  NNEW(eei,double,6*mi[0]**ne);
439  if(*nener==1){
440  NNEW(stiini,double,6*mi[0]**ne);
441  NNEW(enerini,double,mi[0]**ne);}
442 
443  results(co,nk,kon,ipkon,lakon,ne,v,stn,inum,stx,
444  elcon,nelcon,rhcon,nrhcon,alcon,nalcon,alzero,ielmat,
445  ielorien,norien,orab,ntmat_,t0,t1act,ithermal,
446  prestr,iprestr,filab,eme,emn,een,iperturb,
447  f,fn,nactdof,&iout,qa,vold,b,nodeboun,ndirboun,xbounact,nboun,ipompc,
448  nodempc,coefmpc,labmpc,nmpc,nmethod,cam,neq,veold,accold,&bet,
449  &gam,&dtime,&time,ttime,plicon,nplicon,plkcon,nplkcon,
450  xstateini,xstiff,xstate,npmat_,epn,matname,mi,&ielas,&icmd,
451  ncmat_,nstate_,stiini,vini,ikboun,ilboun,ener,enern,emeini,
452  xstaten,eei,enerini,cocon,ncocon,set,nset,istartset,iendset,
453  ialset,nprint,prlab,prset,qfx,qfn,trab,inotr,ntrans,fmpc,
454  nelemload,nload,ikmpc,ilmpc,istep,&iinc,springarea,&reltime,
455  &ne0,xforc,nforc,thicke,shcon,nshcon,
456  sideload,xloadact,xloadold,&icfd,inomat,pslavsurf,pmastsurf,
457  &mortar,islavact,cdn,islavnode,nslavnode,ntie,clearini,
458  islavsurf,ielprop,prop);
459 
460  SFREE(eei);
461  if(*nener==1){
462  SFREE(stiini);SFREE(enerini);}
463 
464  memcpy(&vold[0],&v[0],sizeof(double)*mt**nk);
465  memcpy(&sti[0],&stx[0],sizeof(double)*6*mi[0]**ne);
466 /* for(k=0;k<mt**nk;++k){
467  vold[k]=v[k];
468  }
469  for(k=0;k<6*mi[0]**ne;++k){
470  sti[k]=stx[k];
471  }*/
472 
473  ++*kode;
474 
475  /* for cyclic symmetric sectors: duplicating the results */
476 
477  if(*mcs>0){
478  ptime=*ttime+time;
479  frdcyc(co,nk,kon,ipkon,lakon,ne,v,stn,inum,nmethod,kode,filab,een,t1act,
480  fn,&ptime,epn,ielmat,matname,cs,mcs,nkon,enern,xstaten,
481  nstate_,istep,&iinc,iperturb,ener,mi,output,ithermal,
482  qfn,ialset,istartset,iendset,trab,inotr,ntrans,orab,
483  ielorien,norien,sti,veold,&noddiam,set,nset,emn,thicke,
484  jobnamec,&ne0,cdn,&mortar,nmat);
485  }
486  else{
487  if(strcmp1(&filab[1044],"ZZS")==0){
488  NNEW(neigh,ITG,40**ne);
489  NNEW(ipneigh,ITG,*nk);
490  }
491  ptime=*ttime+time;
492  frd(co,nk,kon,ipkon,lakon,ne,v,stn,inum,nmethod,
493  kode,filab,een,t1act,fn,&ptime,epn,ielmat,matname,enern,xstaten,
494  nstate_,istep,&iinc,ithermal,qfn,&mode,&noddiam,trab,inotr,
495  ntrans,orab,ielorien,norien,description,ipneigh,neigh,
496  mi,stx,vr,vi,stnr,stni,vmax,stnmax,&ngraph,veold,ener,ne,
497  cs,set,nset,istartset,iendset,ialset,eenmax,fnr,fni,emn,
498  thicke,jobnamec,output,qfx,cdn,&mortar,cdnr,cdni,nmat);
499  if(strcmp1(&filab[1044],"ZZS")==0){SFREE(ipneigh);SFREE(neigh);}
500  }
501 
502  SFREE(v);SFREE(stn);SFREE(inum);
503  SFREE(b);SFREE(stx);SFREE(fn);
504 
505  if(strcmp1(&filab[261],"E ")==0) SFREE(een);
506  if(strcmp1(&filab[522],"ENER")==0) SFREE(enern);
507  if(strcmp1(&filab[2697],"ME ")==0) SFREE(emn);
508 
509  }
510  else {
511 
512  /* error occurred in mafill: storing the geometry in frd format */
513 
514  ++*kode;
515  NNEW(inum,ITG,*nk);for(k=0;k<*nk;k++) inum[k]=1;
516  if(strcmp1(&filab[1044],"ZZS")==0){
517  NNEW(neigh,ITG,40**ne);
518  NNEW(ipneigh,ITG,*nk);
519  }
520  ptime=*ttime+time;
521  frd(co,nk,kon,ipkon,lakon,ne,v,stn,inum,nmethod,
522  kode,filab,een,t1,fn,&ptime,epn,ielmat,matname,enern,xstaten,
523  nstate_,istep,&iinc,ithermal,qfn,&mode,&noddiam,trab,inotr,
524  ntrans,orab,ielorien,norien,description,ipneigh,neigh,
525  mi,sti,vr,vi,stnr,stni,vmax,stnmax,&ngraph,veold,ener,ne,
526  cs,set,nset,istartset,iendset,ialset,eenmax,fnr,fni,emn,
527  thicke,jobnamec,output,qfx,cdn,&mortar,cdnr,cdni,nmat);
528  if(strcmp1(&filab[1044],"ZZS")==0){SFREE(ipneigh);SFREE(neigh);}
529  SFREE(inum);FORTRAN(stop,());
530 
531  }
532 
533  /* updating the loading at the end of the step;
534  important in case the amplitude at the end of the step
535  is not equal to one */
536 
537  for(k=0;k<*nboun;++k){xbounold[k]=xbounact[k];}
538  for(k=0;k<*nforc;++k){xforcold[k]=xforcact[k];}
539  for(k=0;k<2**nload;++k){xloadold[k]=xloadact[k];}
540  for(k=0;k<7**nbody;k=k+7){xbodyold[k]=xbodyact[k];}
541  if(*ithermal==1){
542  for(k=0;k<*nk;++k){t1old[k]=t1act[k];}
543  for(k=0;k<*nk;++k){vold[mt*k]=t1act[k];}
544  }
545 
546  SFREE(xbounact);SFREE(xforcact);SFREE(xloadact);SFREE(t1act);SFREE(ampli);
547  SFREE(xbodyact);if(*nbody>0) SFREE(ipobody);SFREE(xstiff);
548 
549  if(iglob==1){SFREE(integerglob);SFREE(doubleglob);}
550 
551  *icolp=icol;
552  *irowp=irow;
553 
554  (*ttime)+=(*tper);
555 
556  return;
557 }
void spooles_solve(double *b, ITG *neq)
void pardiso_main(double *ad, double *au, double *adb, double *aub, double *sigma, double *b, ITG *icol, ITG *irow, ITG *neq, ITG *nzs, ITG *symmetryflag, ITG *inputformat, ITG *jq, ITG *nzs3)
subroutine writesubmatrix(submatrix, noderetain, ndirretain, nretain, jobnamec)
Definition: writesubmatrix.f:19
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 preiter(double *ad, double **aup, double *b, ITG **icolp, ITG **irowp, ITG *neq, ITG *nzs, ITG *isolver, ITG *iperturb)
Definition: preiter.c:23
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))
ITG strcmp1(const char *s1, const char *s2)
Definition: strcmp1.c:24
void frdcyc(double *co, ITG *nk, ITG *kon, ITG *ipkon, char *lakon, ITG *ne, 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 *cs, ITG *mcs, ITG *nkon, double *enern, double *xstaten, ITG *nstate_, ITG *istep, ITG *iinc, ITG *iperturb, double *ener, ITG *mi, char *output, ITG *ithermal, double *qfn, ITG *ialset, ITG *istartset, ITG *iendset, double *trab, ITG *inotr, ITG *ntrans, double *orab, ITG *ielorien, ITG *norien, double *sti, double *veold, ITG *noddiam, char *set, ITG *nset, double *emn, double *thicke, char *jobnamec, ITG *ne0, double *cdn, ITG *mortar, ITG *nmat)
Definition: frdcyc.c:24
#define DMEMSET(a, b, c, d)
Definition: CalculiX.h:45
void spooles(double *ad, double *au, double *adb, double *aub, double *sigma, double *b, ITG *icol, ITG *irow, ITG *neq, ITG *nzs, ITG *symmtryflag, ITG *inputformat, ITG *nzs3)
subroutine stop()
Definition: stop.f:19
void tau(double *ad, double **aup, double *adb, double *aubp, double *sigma, double *b, ITG *icol, ITG **irowp, ITG *neq, ITG *nzs)
void pardiso_factor(double *ad, double *au, double *adb, double *aub, double *sigma, ITG *icol, ITG *irow, ITG *neq, ITG *nzs, ITG *symmetryflag, ITG *inputformat, ITG *jq, ITG *nzs3)
void sgi_main(double *ad, double *au, double *adb, double *aub, double *sigma, double *b, ITG *icol, ITG *irow, ITG *neq, ITG *nzs, ITG token)
#define RENEW(a, b, c)
Definition: CalculiX.h:40
#define SFREE(a)
Definition: CalculiX.h:41
void results(double *co, ITG *nk, ITG *kon, ITG *ipkon, char *lakon, ITG *ne, double *v, double *stn, ITG *inum, double *stx, double *elcon, ITG *nelcon, double *rhcon, ITG *nrhcon, double *alcon, ITG *nalcon, double *alzero, ITG *ielmat, ITG *ielorien, ITG *norien, double *orab, ITG *ntmat_, double *t0, double *t1, ITG *ithermal, double *prestr, ITG *iprestr, char *filab, double *eme, double *emn, double *een, ITG *iperturb, double *f, double *fn, ITG *nactdof, ITG *iout, double *qa, double *vold, double *b, ITG *nodeboun, ITG *ndirboun, double *xboun, ITG *nboun, ITG *ipompc, ITG *nodempc, double *coefmpc, char *labmpc, ITG *nmpc, ITG *nmethod, double *vmax, ITG *neq, double *veold, double *accold, double *beta, double *gamma, double *dtime, double *time, double *ttime, double *plicon, ITG *nplicon, double *plkcon, ITG *nplkcon, double *xstateini, double *xstiff, double *xstate, ITG *npmat_, double *epl, char *matname, ITG *mi, ITG *ielas, ITG *icmd, ITG *ncmat_, ITG *nstate_, double *stiini, double *vini, ITG *ikboun, ITG *ilboun, double *ener, double *enern, double *emeini, double *xstaten, double *eei, double *enerini, double *cocon, ITG *ncocon, char *set, ITG *nset, ITG *istartset, ITG *iendset, ITG *ialset, ITG *nprint, char *prlab, char *prset, double *qfx, double *qfn, double *trab, ITG *inotr, ITG *ntrans, double *fmpc, ITG *nelemload, ITG *nload, ITG *ikmpc, ITG *ilmpc, ITG *istep, ITG *iinc, double *springarea, double *reltime, ITG *ne0, double *xforc, ITG *nforc, double *thicke, double *shcon, ITG *nshcon, char *sideload, double *xload, double *xloadold, ITG *icfd, ITG *inomat, double *pslavsurf, double *pmastsurf, ITG *mortar, ITG *islavact, double *cdn, ITG *islavnode, ITG *nslavnode, ITG *ntie, double *clearini, ITG *islavsurf, ITG *ielprop, double *prop)
Definition: results.c:42
subroutine tempload(xforcold, xforc, xforcact, iamforc, nforc, xloadold, xload, xloadact, iamload, nload, ibody, xbody, nbody, xbodyold, xbodyact, t1old, t1, t1act, iamt1, nk, amta, namta, nam, ampli, time, reltime, ttime, dtime, ithermal, nmethod, xbounold, xboun, xbounact, iamboun, nboun, nodeboun, ndirboun, nodeforc, ndirforc, istep, iinc, co, vold, itg, ntg, amname, ikboun, ilboun, nelemload, sideload, mi, ntrans, trab, inotr, veold, integerglob, doubleglob, tieset, istartset, iendset, ialset, ntie, nmpc, ipompc, ikmpc, ilmpc, nodempc, coefmpc)
Definition: tempload.f:19
subroutine mafillsm(co, nk, kon, ipkon, lakon, ne, nodeboun, ndirboun, xboun, nboun, ipompc, nodempc, coefmpc, nmpc, nodeforc, ndirforc, xforc, nforc, nelemload, sideload, xload, nload, xbody, ipobody, nbody, cgr, ad, au, fext, nactdof, icol, jq, irow, neq, nzl, nmethod, ikmpc, ilmpc, ikboun, ilboun, elcon, nelcon, rhcon, nrhcon, alcon, nalcon, alzero, ielmat, ielorien, norien, orab, ntmat_, t0, t1, ithermal, prestr, iprestr, vold, iperturb, sti, nzs, stx, adb, aub, iexpl, plicon, nplicon, plkcon, nplkcon, xstiff, npmat_, dtime, matname, mi, ncmat_, mass, stiffness, buckling, rhsi, intscheme, physcon, shcon, nshcon, cocon, ncocon, ttime, time, istep, iinc, coriolis, ibody, xloadold, reltime, veold, springarea, nstate_, xstateini, xstate, thicke, integerglob, doubleglob, tieset, istartset, iendset, ialset, ntie, nasym, pslavsurf, pmastsurf, mortar, clearini, ielprop, prop)
Definition: mafillsm.f:19
void getglobalresults(char *jobnamec, ITG **integerglobp, double **doubleglobp, ITG *nboun, ITG *iamboun, double *xboun, ITG *nload, char *sideload, ITG *iamload, ITG *iglob, ITG *nforc, ITG *iamforc, double *xforc, ITG *ithermal, ITG *nk, double *t1, ITG *iamt1)
void spooles_factor(double *ad, double *au, double *adb, double *aub, double *sigma, ITG *icol, ITG *irow, ITG *neq, ITG *nzs, ITG *symmetryflag, ITG *inputformat, ITG *nzs3)
#define ITG
Definition: CalculiX.h:51
void pardiso_solve(double *b, ITG *neq, ITG *symmetryflag)
#define NNEW(a, b, c)
Definition: CalculiX.h:39
subroutine bodyforce(cbody, ibody, ipobody, nbody, set, istartset, iendset, ialset, inewton, nset, ifreebody, k)
Definition: bodyforce.f:19
Hosted by OpenAircraft.com, (Michigan UAV, LLC)