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

Go to the source code of this file.

Functions

int main (int argc, char *argv[])
 

Function Documentation

int main ( int  argc,
char *  argv[] 
)

Definition at line 36 of file ccx_2.8.c.

37 {
38 
39 FILE *f1;
40 
41 ITG *kon=NULL, *nodeboun=NULL, *ndirboun=NULL, *ipompc=NULL,
42  *nodempc=NULL, *nodeforc=NULL, *ndirforc=NULL,
43  *nelemload=NULL,im,*inodesd=NULL,nload1,*idefforc=NULL,
44  *nactdof=NULL, *icol=NULL,*ics=NULL,
45  *jq=NULL, *mast1=NULL, *irow=NULL, *rig=NULL,*idefbody=NULL,
46  *ikmpc=NULL, *ilmpc=NULL, *ikboun=NULL, *ilboun=NULL,
47  *nreorder=NULL,*ipointer=NULL,*idefload=NULL,
48  *istartset=NULL, *iendset=NULL, *ialset=NULL, *ielmat=NULL,
49  *ielorien=NULL, *nrhcon=NULL, *nodebounold=NULL, *ndirbounold=NULL,
50  *nelcon=NULL, *nalcon=NULL, *iamforc=NULL, *iamload=NULL,
51  *iamt1=NULL, *namta=NULL, *ipkon=NULL, *iamboun=NULL,
52  *nplicon=NULL, *nplkcon=NULL, *inotr=NULL, *iponor=NULL, *knor=NULL,
53  *ikforc=NULL, *ilforc=NULL, *iponoel=NULL, *inoel=NULL, *nshcon=NULL,
54  *ncocon=NULL,*ibody=NULL,*ielprop=NULL,*islavsurf=NULL,
55  *ipoinpc=NULL,icfd=0,mt,nxstate,nload0,iload;
56 
57 double *co=NULL, *xboun=NULL, *coefmpc=NULL, *xforc=NULL,*clearini=NULL,
58  *xload=NULL, *xbounold=NULL, *xforcold=NULL,
59  *vold=NULL, *sti=NULL, *xloadold=NULL, *xnor=NULL,
60  *reorder=NULL,*dcs=NULL, *thickn=NULL, *thicke=NULL, *offset=NULL,
61  *elcon=NULL, *rhcon=NULL, *alcon=NULL, *alzero=NULL, *t0=NULL, *t1=NULL,
62  *prestr=NULL, *orab=NULL, *amta=NULL, *veold=NULL, *accold=NULL,
63  *t1old=NULL, *eme=NULL, *plicon=NULL, *pslavsurf=NULL, *plkcon=NULL,
64  *xstate=NULL, *trab=NULL, *ener=NULL, *shcon=NULL, *cocon=NULL,
65  *cs=NULL,*tietol=NULL,*fmpc=NULL,*prop=NULL,*t0g=NULL,*t1g=NULL,
66  *xbody=NULL,*xbodyold=NULL;
67 
68 double ctrl[27]={4.5,8.5,9.5,16.5,10.5,4.5,0.,5.5,0.,0.,0.25,0.5,0.75,0.85,0.,0.,1.5,0.,0.005,0.01,0.,0.,0.02,1.e-5,1.e-3,1.e-8,1.e30};
69 
70 char *sideload=NULL, *set=NULL, *matname=NULL, *orname=NULL, *amname=NULL,
71  *filab=NULL, *lakon=NULL, *labmpc=NULL, *prlab=NULL, *prset=NULL,
72  jobnamec[660]="",jobnamef[132]="",output[4]="asc", *typeboun=NULL,
73  *inpc=NULL,*tieset=NULL,*cbody=NULL,fneig[132]="",*sideloadtemp=NULL,
74  kind1[2]="T",kind2[2]="T",*heading=NULL;
75 
76 ITG nk,ne,nboun,nmpc,nforc,nload,nprint,nset,nalset,nentries=15,
77  nmethod,neq[3]={0,0,0},i,mpcfree=1,mei[4],j,nzl,nam,nbounold=0,
78  nforcold=0,nloadold=0,nbody,nbody_=0,nbodyold=0,network=0,nheading_=0,
79  k,nzs[3],nmpc_=0,nload_=0,nforc_=0,istep,istat,nboun_=0,nintpoint=0,
80  iperturb[2]={0,0},nmat,ntmat_=0,norien,ithermal[2]={0,0},nmpcold,
81  iprestr,kode,isolver=0,nslavs=0,nkon_=0,ne0,nkon0,mortar=0,
82  jout[2]={1,1},nlabel,nkon=0,idrct,jmax[2],iexpl,nevtot=0,ifacecount=0,
83  iplas=0,npmat_=0,mi[3]={0,3,1},ntrans,mpcend=-1,namtot_=0,iumat=0,mpcmult,
84  icascade=0,maxlenmpc,mpcinfo[4],ne1d=0,ne2d=0,infree[4]={0,0,0,0},
85  callfrommain,nflow=0,jin=0,irstrt=0,nener=0,jrstrt=0,nenerold,
86  nline,ipoinp[2*nentries],*inp=NULL,ntie,ntie_=0,mcs=0,nprop_=0,
87  nprop=0,itpamp=0,iviewfile,nkold,nevdamp_=0,npt_=0,cyclicsymmetry,
88  nmethodl,iaxial=1,inext;
89 
90 ITG *meminset=NULL,*rmeminset=NULL;
91 
92 ITG nzs_,nk_=0,ne_=0,nset_=0,nalset_=0,nmat_=0,norien_=0,nam_=0,
93  ntrans_=0,ncs_=0,nstate_=0,ncmat_=0,memmpc_=0,nprint_=0;
94 
95 double fei[3],*xmodal=NULL,timepar[5],
96  alpha,ttime=0.,qaold[2]={0.,0.},physcon[10]={0.,0.,0.,0.,0.,0.,0.,0.,0.,0.};
97 
98 #ifdef CALCULIX_MPI
99 MPI_Init(&argc, &argv) ;
100 MPI_Comm_rank(MPI_COMM_WORLD, &myid) ;
101 MPI_Comm_size(MPI_COMM_WORLD, &nproc) ;
102 #endif
103 
104 if(argc==1){printf("Usage: CalculiX.exe -i jobname\n");FORTRAN(stop,());}
105 else{
106  for(i=1;i<argc;i++){
107  if(strcmp1(argv[i],"-i")==0) {
108  strcpy(jobnamec,argv[i+1]);strcpy1(jobnamef,argv[i+1],132);jin++;break;}
109  if(strcmp1(argv[i],"-v")==0) {
110  printf("\nThis is version version 2.8p1\n\n");
111  FORTRAN(stop,());
112  }
113  }
114  if(jin==0){strcpy(jobnamec,argv[1]);strcpy1(jobnamef,argv[1],132);}
115 
116  for(i=1;i<argc;i++){
117  if(strcmp1(argv[i],"-o")==0) {
118  strcpy(output,argv[i+1]);break;}
119  }
120 }
121 
122 #ifdef BAM
123 ITG lop=0,lrestart=0,kstep=1,kinc=1;
124 double time[2],dtime;
125 FORTRAN(uexternaldb,(&lop,&lrestart,time,&dtime,&kstep,&kinc));
126 #endif
127 
128 FORTRAN(openfile,(jobnamef,output));
129 
130 printf("\n************************************************************\n\n");
131 printf("CalculiX version 2.8p1, Copyright(C) 1998-2015 Guido Dhondt\n");
132 printf("CalculiX Extras version 2.8p1, Copyright(C) 2013-2015 Peter Gustafson\n");
133 printf("CalculiX comes with ABSOLUTELY NO WARRANTY. This is free\n");
134 printf("software, and you are welcome to redistribute it under\n");
135 printf("certain conditions, see gpl.htm\n\n");
136 printf("************************************************************\n\n");
137 printf("You are using an executable made on So 8. Feb 13:14:01 CET 2015\n");
138 fflush(stdout);
139 
140 istep=0;
141 istat=0;
142 iprestr=0;
143 kode=0;
144 
145 /* default solver */
146 
147 #if defined(SGI)
148  isolver=4;
149 #elif defined(PARDISO)
150  isolver=7;
151 #elif defined(SPOOLES)
152  isolver=0;
153 #elif defined(TAUCS)
154  isolver=5;
155 #else
156  isolver=3;
157 #endif
158 
159 /* conservative estimate of the fields to be allocated */
160 
161 readinput(jobnamec,&inpc,&nline,&nset_,ipoinp,&inp,&ipoinpc,ithermal);
162 
163 NNEW(set,char,81*nset_);
164 NNEW(meminset,ITG,nset_);
165 NNEW(rmeminset,ITG,nset_);
166 
167 FORTRAN(allocation,(&nload_,&nforc_,&nboun_,&nk_,&ne_,&nmpc_,&nset_,&nalset_,
168  &nmat_,&ntmat_,&npmat_,&norien_,&nam_,&nprint_,mi,&ntrans_,
169  set,meminset,rmeminset,&ncs_,&namtot_,&ncmat_,&memmpc_,&ne1d,
170  &ne2d,&nflow,jobnamec,&irstrt,ithermal,&nener,&nstate_,&istep,
171  inpc,ipoinp,inp,&ntie_,&nbody_,&nprop_,ipoinpc,&nevdamp_,&npt_,&nslavs,
172  &nkon_,&mcs,&mortar,&ifacecount,&nintpoint,infree,&nheading_));
173 
174 SFREE(set);SFREE(meminset);SFREE(rmeminset);mt=mi[1]+1;
175 NNEW(heading,char,66*nheading_);
176 
177 nzs_=20000000;
178 
179 nload=0;nbody=0;nforc=0;nboun=0;nk=0;nmpc=0;nam=0;
180 
181 /* caveat: if changing next line:
182  - change noelfiles appropriately
183  - change nlabel in geomview.f, expand.c, storeresidual.f
184  and createmddof.f
185  - change the dimension of label in geomview.f */
186 
187 nlabel=46;
188 
189 while(istat>=0) {
190 
191  fflush(stdout);
192 
193  /* in order to reduce the number of variables to be transferred to
194  the subroutines, the max. field sizes are (for most fields) copied
195  into the real sizes */
196 
197  nzs[1]=nzs_;
198  nprint=nprint_;
199 
200  if((istep==0)||(irstrt<0)) {
201  ne=ne_;
202  nset=nset_;
203  nalset=nalset_;
204  nmat=nmat_;
205  norien=norien_;
206  ntrans=ntrans_;
207  ntie=ntie_;
208 
209  /* allocating space before the first step */
210 
211  /* coordinates and topology */
212 
213  NNEW(co,double,3*nk_);
214  NNEW(kon,ITG,nkon_);
215  NNEW(ipkon,ITG,ne_);
216  NNEW(lakon,char,8*ne_);
217 
218  /* property cards */
219 
220  if(nprop_>0){
221  NNEW(ielprop,ITG,ne_);
222  for(i=0;i<ne_;i++) ielprop[i]=-1;
223  NNEW(prop,double,nprop_);
224  }
225 
226  /* fields for 1-D and 2-D elements */
227 
228  if((ne1d!=0)||(ne2d!=0)){
229  NNEW(iponor,ITG,2*nkon_);
230  for(i=0;i<2*nkon_;i++) iponor[i]=-1;
231  NNEW(xnor,double,36*ne1d+24*ne2d);
232  NNEW(knor,ITG,24*(ne1d+ne2d)*(mi[2]+1));
233  NNEW(thickn,double,2*nk_);
234  NNEW(thicke,double,mi[2]*nkon_);
235  NNEW(offset,double,2*ne_);
236  NNEW(iponoel,ITG,nk_);
237  NNEW(inoel,ITG,9*ne1d+24*ne2d);
238  NNEW(rig,ITG,nk_);
239  if(infree[2]==0)infree[2]=1;
240  }
241 
242  /* SPC's */
243 
244  NNEW(nodeboun,ITG,nboun_);
245  NNEW(ndirboun,ITG,nboun_);
246  NNEW(typeboun,char,nboun_+1);
247  if((istep == 0)||((irstrt<0)&&(nam_>0)))NNEW(iamboun,ITG,nboun_);
248  NNEW(xboun,double,nboun_);
249  NNEW(ikboun,ITG,nboun_);
250  NNEW(ilboun,ITG,nboun_);
251 
252  /* MPC's */
253 
254  NNEW(ipompc,ITG,nmpc_);
255  NNEW(nodempc,ITG,3*memmpc_);
256  for(i=0;i<3*memmpc_;i+=3){nodempc[i+2]=i/3+2;}
257  nodempc[3*memmpc_-1]=0;
258  NNEW(coefmpc,double,memmpc_);
259  NNEW(labmpc,char,20*nmpc_+1);
260  NNEW(ikmpc,ITG,nmpc_);
261  NNEW(ilmpc,ITG,nmpc_);
262  NNEW(fmpc,double,nmpc_);
263 
264  /* nodal loads */
265 
266  NNEW(nodeforc,ITG,2*nforc_);
267  NNEW(ndirforc,ITG,nforc_);
268  if((istep == 0)||((irstrt<0)&&(nam_>0)))NNEW(iamforc,ITG,nforc_);
269  NNEW(idefforc,ITG,nforc_);
270  NNEW(xforc,double,nforc_);
271  NNEW(ikforc,ITG,nforc_);
272  NNEW(ilforc,ITG,nforc_);
273 
274  /* distributed facial loads */
275 
276  NNEW(nelemload,ITG,2*nload_);
277  if((istep == 0)||((irstrt<0)&&(nam_>0)))NNEW(iamload,ITG,2*nload_);
278  NNEW(idefload,ITG,nload_);
279  NNEW(sideload,char,20*nload_);
280  NNEW(xload,double,2*nload_);
281 
282  /* distributed volumetric loads */
283 
284  NNEW(cbody,char,81*nbody_);
285  NNEW(idefbody,ITG,nbody_);
286  NNEW(ibody,ITG,3*nbody_);
287  NNEW(xbody,double,7*nbody_);
288  NNEW(xbodyold,double,7*nbody_);
289 
290  /* printing output */
291 
292  NNEW(prlab,char,6*nprint_);
293  NNEW(prset,char,81*nprint_);
294 
295  /* set definitions */
296 
297  NNEW(set,char,81*nset);
298  NNEW(istartset,ITG,nset);
299  NNEW(iendset,ITG,nset);
300  NNEW(ialset,ITG,nalset);
301 
302  /* (hyper)elastic constants */
303 
304  NNEW(elcon,double,(ncmat_+1)*ntmat_*nmat);
305  NNEW(nelcon,ITG,2*nmat);
306 
307  /* density */
308 
309  NNEW(rhcon,double,2*ntmat_*nmat);
310  NNEW(nrhcon,ITG,nmat);
311 
312  /* specific heat */
313 
314  NNEW(shcon,double,4*ntmat_*nmat);
315  NNEW(nshcon,ITG,nmat);
316 
317  /* thermal expansion coefficients */
318 
319  NNEW(alcon,double,7*ntmat_*nmat);
320  NNEW(nalcon,ITG,2*nmat);
321  NNEW(alzero,double,nmat);
322 
323  /* conductivity */
324 
325  NNEW(cocon,double,7*ntmat_*nmat);
326  NNEW(ncocon,ITG,2*nmat);
327 
328  /* isotropic and kinematic hardening coefficients*/
329 
330  if(npmat_>0){
331  NNEW(plicon,double,(2*npmat_+1)*ntmat_*nmat);
332  NNEW(nplicon,ITG,(ntmat_+1)*nmat);
333  NNEW(plkcon,double,(2*npmat_+1)*ntmat_*nmat);
334  NNEW(nplkcon,ITG,(ntmat_+1)*nmat);
335  }
336 
337  /* linear dynamic properties */
338 
339  NNEW(xmodal,double,11+nevdamp_);
340  xmodal[10]=nevdamp_+0.5;
341 
342  /* internal state variables (nslavs is needed for restart
343  calculations) */
344 
345  if(mortar!=1){
346  NNEW(xstate,double,nstate_*mi[0]*(ne+nslavs));
347  nxstate=nstate_*mi[0]*(ne+nslavs);
348  }else if(mortar==1){
349  NNEW(xstate,double,nstate_*mi[0]*(ne+nintpoint));
350  nxstate=nstate_*mi[0]*(ne+nintpoint);
351  }
352 
353  /* material orientation */
354 
355  if((istep == 0)||((irstrt<0)&&(norien>0))) {
356  NNEW(orname,char,80*norien);
357  NNEW(orab,double,7*norien);
358  NNEW(ielorien,ITG,mi[2]*ne_);
359  }
360 
361  /* transformations */
362 
363  if((istep == 0)||((irstrt<0)&&(ntrans>0))) {
364  NNEW(trab,double,7*ntrans);
365  NNEW(inotr,ITG,2*nk_);
366  }
367 
368  /* amplitude definitions */
369 
370  if((istep == 0)||((irstrt<0)&&(nam_>0))) {
371  NNEW(amname,char,80*nam_);
372  NNEW(amta,double,2*namtot_);
373  NNEW(namta,ITG,3*nam_);
374  }
375 
376  if((istep == 0)||((irstrt<0)&&(ithermal[0]>0))) {
377  NNEW(t0,double,nk_);
378  NNEW(t1,double,nk_);
379  if((ne1d!=0)||(ne2d!=0)){
380  NNEW(t0g,double,2*nk_);
381  NNEW(t1g,double,2*nk_);
382  }
383  }
384 
385  /* the number in next line is NOT 1.2357111317 -> points
386  to user input; instead it is a generic nonzero
387  initialization */
388 
389  if(istep==0){
390  DMEMSET(t0,0,nk_,1.2357111319);
391  DMEMSET(t1,0,nk_,1.2357111319);
392  }
393 
394  if((istep == 0)||((irstrt<0)&&(ithermal[0]>0)&&(nam_>0)))NNEW(iamt1,ITG,nk_);
395 
396  if((istep==0)||((irstrt<0)&&(iprestr>0)))NNEW(prestr,double,6*mi[0]*ne_);
397 
398  NNEW(vold,double,mt*nk_);
399  NNEW(veold,double,mt*nk_);
400 
401  NNEW(ielmat,ITG,mi[2]*ne_);
402 
403  NNEW(matname,char,80*nmat);
404 
405  NNEW(filab,char,87*nlabel);
406 
407  /* tied constraints */
408 
409  if(ntie_>0){
410  NNEW(tieset,char,243*ntie_);
411  NNEW(tietol,double,3*ntie_);
412  NNEW(cs,double,17*ntie_);
413  }
414 
415  /* temporary fields for cyclic symmetry calculations */
416 
417  if((ncs_>0)||(npt_>0)){
418  if(2*npt_>24*ncs_){
419  NNEW(ics,ITG,2*npt_);
420  }else{
421  NNEW(ics,ITG,24*ncs_);
422  }
423  if(npt_>30*ncs_){
424  NNEW(dcs,double,npt_);
425  }else{
426  NNEW(dcs,double,30*ncs_);
427  }
428  }
429 
430  /* slave faces */
431 
432  NNEW(islavsurf,ITG,2*ifacecount+2);
433 
434  }
435  else {
436 
437  /* allocating and reallocating space for subsequent steps */
438 
439  if((nmethod!=4)&&(nmethod!=5)&&(nmethod!=8)&&(nmethod!=9)&&
440  ((abs(nmethod)!=1)||(iperturb[0]<2))){
441  NNEW(veold,double,mt*nk_);
442  }
443  else{
444  RENEW(veold,double,mt*nk_);
445  DMEMSET(veold,mt*nk,mt*nk_,0.);
446  }
447  RENEW(vold,double,mt*nk_);
448  DMEMSET(vold,mt*nk,mt*nk_,0.);
449 
450  RENEW(nodeboun,ITG,nboun_);
451  RENEW(ndirboun,ITG,nboun_);
452  RENEW(typeboun,char,nboun_+1);
453  RENEW(xboun,double,nboun_);
454  RENEW(ikboun,ITG,nboun_);
455  RENEW(ilboun,ITG,nboun_);
456 
457  RENEW(nodeforc,ITG,2*nforc_);
458  RENEW(ndirforc,ITG,nforc_);
459  NNEW(idefforc,ITG,nforc_);
460  RENEW(xforc,double,nforc_);
461  RENEW(ikforc,ITG,nforc_);
462  RENEW(ilforc,ITG,nforc_);
463 
464  RENEW(nelemload,ITG,2*nload_);
465  NNEW(idefload,ITG,nload_);
466  RENEW(sideload,char,20*nload_);
467  RENEW(xload,double,2*nload_);
468 
469  RENEW(cbody,char,81*nbody_);
470  NNEW(idefbody,ITG,nbody_);
471  RENEW(ibody,ITG,3*nbody_);
472  RENEW(xbody,double,7*nbody_);
473  RENEW(xbodyold,double,7*nbody_);
474  for(i=7*nbodyold;i<7*nbody_;i++) xbodyold[i]=0;
475 
476  if(nam > 0) {
477  RENEW(iamforc,ITG,nforc_);
478  RENEW(iamload,ITG,2*nload_);
479  RENEW(iamboun,ITG,nboun_);
480  RENEW(amname,char,80*nam_);
481  RENEW(amta,double,2*namtot_);
482  RENEW(namta,ITG,3*nam_);
483  }
484 
485  RENEW(ipompc,ITG,nmpc_);
486 
487  RENEW(labmpc,char,20*nmpc_+1);
488  RENEW(ikmpc,ITG,nmpc_);
489  RENEW(ilmpc,ITG,nmpc_);
490  RENEW(fmpc,double,nmpc_);
491 
492  if(ntrans > 0){
493  RENEW(inotr,ITG,2*nk_);DMEMSET(inotr,2*nk,2*nk_,0);
494  }
495 
496  RENEW(co,double,3*nk_);DMEMSET(co,3*nk,3*nk_,0.);
497 
498  if(ithermal[0] != 0){
499  RENEW(t0,double,nk_);DMEMSET(t0,nk,nk_,0.);
500  RENEW(t1,double,nk_);DMEMSET(t1,nk,nk_,0.);
501  if((ne1d!=0)||(ne2d!=0)){
502  RENEW(t0g,double,2*nk_);DMEMSET(t0g,2*nk,2*nk_,0.);
503  RENEW(t1g,double,2*nk_);DMEMSET(t1g,2*nk,2*nk_,0.);
504  }
505  if(nam > 0) {RENEW(iamt1,ITG,nk_);}
506  }
507 
508  }
509 
510  /* allocation of fields in the restart file */
511 
512  if(irstrt<0){
513  NNEW(nodebounold,ITG,nboun_);
514  NNEW(ndirbounold,ITG,nboun_);
515  NNEW(xbounold,double,nboun_);
516  NNEW(xforcold,double,nforc_);
517  NNEW(xloadold,double,2*nload_);
518  if(ithermal[0]!=0) NNEW(t1old,double,nk_);
519  NNEW(sti,double,6*mi[0]*ne);
520  NNEW(eme,double,6*mi[0]*ne);
521  if(nener==1)NNEW(ener,double,mi[0]*ne*2);
522  if(mcs>ntie_) RENEW(cs,double,17*mcs);
523  if(mortar==1){
524  NNEW(pslavsurf,double,3*nintpoint);
525  NNEW(clearini,double,3*9*ifacecount);
526  }
527 
528  }
529 
530  nenerold=nener;
531  nkold=nk;
532 
533  /* opening the eigenvalue file and checking for cyclic symmetry */
534 
535  strcpy(fneig,jobnamec);
536  strcat(fneig,".eig");
537  cyclicsymmetry=0;
538  if((f1=fopen(fneig,"rb"))!=NULL){
539  if(fread(&cyclicsymmetry,sizeof(ITG),1,f1)!=1){
540  printf("*ERROR reading the information whether cyclic symmetry is involved in the eigenvalue file");
541  exit(0);
542  }
543  fclose(f1);
544  }
545 
546  nmpcold=nmpc;
547 
548  /* reading the input file */
549 
550  FORTRAN(calinput,(co,&nk,kon,ipkon,lakon,&nkon,&ne,
551  nodeboun,ndirboun,xboun,&nboun,
552  ipompc,nodempc,coefmpc,&nmpc,&nmpc_,nodeforc,ndirforc,xforc,&nforc,
553  &nforc_,nelemload,sideload,xload,&nload,&nload_,
554  &nprint,prlab,prset,&mpcfree,&nboun_,mei,set,istartset,iendset,
555  ialset,&nset,&nalset,elcon,nelcon,rhcon,nrhcon,alcon,nalcon,
556  alzero,t0,t1,matname,ielmat,orname,orab,ielorien,amname,
557  amta,namta,&nam,&nmethod,iamforc,iamload,iamt1,
558  ithermal,iperturb,&istat,&istep,&nmat,&ntmat_,&norien,prestr,
559  &iprestr,&isolver,fei,veold,timepar,
560  xmodal,filab,jout,&nlabel,&idrct,
561  jmax,&iexpl,&alpha,iamboun,plicon,nplicon,
562  plkcon,nplkcon,&iplas,&npmat_,mi,&nk_,trab,inotr,&ntrans,
563  ikboun,ilboun,ikmpc,ilmpc,ics,dcs,&ncs_,&namtot_,cs,&nstate_,
564  &ncmat_,&iumat,&mcs,labmpc,iponor,xnor,knor,thickn,thicke,
565  ikforc,ilforc,offset,iponoel,inoel,rig,infree,nshcon,shcon,
566  cocon,ncocon,physcon,&nflow,
567  ctrl,&maxlenmpc,&ne1d,&ne2d,&nener,vold,nodebounold,
568  ndirbounold,xbounold,xforcold,xloadold,t1old,eme,
569  sti,ener,xstate,jobnamec,&irstrt,&ttime,
570  qaold,output,typeboun,inpc,ipoinp,inp,tieset,tietol,
571  &ntie,fmpc,cbody,ibody,xbody,&nbody,&nbody_,xbodyold,&nam_,
572  ielprop,&nprop,&nprop_,prop,&itpamp,&iviewfile,ipoinpc,&icfd,
573  &nslavs,t0g,t1g,&network,&cyclicsymmetry,idefforc,idefload,
574  idefbody,&mortar,&ifacecount,islavsurf,pslavsurf,clearini,
575  heading,&iaxial));
576 
577  nload0=nload;SFREE(idefforc);SFREE(idefload);SFREE(idefbody);
578 
579  if(nheading_>=0){
580  writeheading(jobnamec,heading,&nheading_);
581  SFREE(heading);
582  nheading_=-1;
583  }
584 
585  if((abs(nmethod)!=1)||(iperturb[0]<2))icascade=0;
586 
587 // FORTRAN(writeboun,(nodeboun,ndirboun,xboun,typeboun,&nboun));
588 
589  if(istat<0) break;
590 
591  if(istep == 1) {
592 
593  /* tied contact constraints: generate appropriate MPC's */
594 
595  tiedcontact(&ntie, tieset, &nset, set,istartset, iendset, ialset,
596  lakon, ipkon, kon,tietol,&nmpc, &mpcfree, &memmpc_,
597  &ipompc, &labmpc, &ikmpc, &ilmpc,&fmpc, &nodempc, &coefmpc,
598  ithermal, co, vold,&icfd,&nmpc_,mi,&nk,&istep,ikboun,&nboun,
599  kind1,kind2);
600 
601  /* reallocating space in the first step */
602 
603  /* allocating and initializing fields pointing to the previous step */
604 
605  RENEW(vold,double,mt*nk);
606  NNEW(sti,double,6*mi[0]*ne);
607 
608  /* strains */
609 
610  NNEW(eme,double,6*mi[0]*ne);
611 
612  /* residual stresses/strains */
613 
614  if(iprestr==1) {
615  RENEW(prestr,double,6*mi[0]*ne);
616  for(i=0;i<ne;i++){
617  for(j=0;j<mi[0];j++){
618  for(k=0;k<6;k++){
619  sti[6*mi[0]*i+6*j+k]=prestr[6*mi[0]*i+6*j+k];
620  }
621  }
622  }
623  }
624  else if(iprestr==2){
625  RENEW(prestr,double,6*mi[0]*ne);
626  for(i=0;i<ne;i++){
627  for(j=0;j<mi[0];j++){
628  for(k=0;k<6;k++){
629  eme[6*mi[0]*i+6*j+k]=prestr[6*mi[0]*i+6*j+k];
630  }
631  }
632  }
633  }
634  else {
635  SFREE(prestr);
636  }
637 
638  NNEW(nodebounold,ITG,nboun);
639  NNEW(ndirbounold,ITG,nboun);
640  NNEW(xbounold,double,nboun);
641  NNEW(xforcold,double,nforc);
642  NNEW(xloadold,double,2*nload);
643 
644  /* initial temperatures: store in the "old" boundary conditions */
645 
646  if(ithermal[0]>1){
647  for(i=0;i<nboun;i++){
648  if(strcmp1(&typeboun[i],"F")==0) continue;
649  if(ndirboun[i]==0){
650  xbounold[i]=vold[mt*(nodeboun[i]-1)];
651  }
652  }
653  }
654 
655  /* initial temperatures: store in the "old" temperature field */
656 
657  if(ithermal[0]!=0){
658  NNEW(t1old,double,nk);
659  for(i=0;i<nk;i++) t1old[i]=t0[i];
660  }
661 
662  /* element definition */
663 
664  RENEW(kon,ITG,nkon);
665  RENEW(ipkon,ITG,ne);
666  RENEW(lakon,char,8*ne);
667 
668  /* property cards */
669 
670  if(nprop>0){
671  RENEW(ielprop,ITG,ne);
672  RENEW(prop,double,nprop);
673  }else{
674  SFREE(ielprop);SFREE(prop);
675  }
676 
677  /* fields for 1-D and 2-D elements */
678 
679  if((ne1d!=0)||(ne2d!=0)){
680  RENEW(iponor,ITG,2*nkon);
681  RENEW(xnor,double,infree[0]);
682  RENEW(knor,ITG,infree[1]);
683  SFREE(thickn);
684  RENEW(thicke,double,mi[2]*nkon);
685  RENEW(offset,double,2*ne);
686  RENEW(inoel,ITG,3*(infree[2]-1));
687  RENEW(iponoel,ITG,infree[3]);
688  RENEW(rig,ITG,infree[3]);
689  }
690 
691  /* set definitions */
692 
693  RENEW(set,char,81*nset);
694  RENEW(istartset,ITG,nset);
695  RENEW(iendset,ITG,nset);
696  RENEW(ialset,ITG,nalset);
697 
698  /* material properties */
699 
700  RENEW(elcon,double,(ncmat_+1)*ntmat_*nmat);
701  RENEW(nelcon,ITG,2*nmat);
702 
703  RENEW(rhcon,double,2*ntmat_*nmat);
704  RENEW(nrhcon,ITG,nmat);
705 
706  RENEW(shcon,double,4*ntmat_*nmat);
707  RENEW(nshcon,ITG,nmat);
708 
709  RENEW(cocon,double,7*ntmat_*nmat);
710  RENEW(ncocon,ITG,2*nmat);
711 
712  RENEW(alcon,double,7*ntmat_*nmat);
713  RENEW(nalcon,ITG,2*nmat);
714  RENEW(alzero,double,nmat);
715 
716  RENEW(matname,char,80*nmat);
717  RENEW(ielmat,ITG,mi[2]*ne);
718 
719  /* allocating space for the state variables */
720 
721  if(mortar!=1){
722  RENEW(xstate,double,nstate_*mi[0]*(ne+nslavs));
723  for(i=nxstate;i<nstate_*mi[0]*(ne+nslavs);i++){xstate[i]=0.;}
724  }else if(mortar==1){
725  RENEW(xstate,double,nstate_*mi[0]*(ne+nintpoint));
726  for(i=nxstate;i<nstate_*mi[0]*(ne+nintpoint);i++){xstate[i]=0.;}
727  }
728 
729  /* next statements for plastic materials and nonlinear springs */
730 
731  if(npmat_>0){
732  RENEW(plicon,double,(2*npmat_+1)*ntmat_*nmat);
733  RENEW(nplicon,ITG,(ntmat_+1)*nmat);
734  RENEW(plkcon,double,(2*npmat_+1)*ntmat_*nmat);
735  RENEW(nplkcon,ITG,(ntmat_+1)*nmat);
736  }
737 
738  /* material orientation */
739 
740  if(norien > 0) {
741  RENEW(orname,char,80*norien);
742  RENEW(ielorien,ITG,mi[2]*ne);
743  RENEW(orab,double,7*norien);
744  }
745  else {
746  SFREE(orname);
747  SFREE(ielorien);
748  SFREE(orab);
749  }
750 
751  /* amplitude definitions */
752 
753  if(nam > 0) {
754  RENEW(amname,char,80*nam);
755  RENEW(namta,ITG,3*nam);
756  RENEW(amta,double,2*namta[3*nam-2]);
757  }
758  else {
759  SFREE(amname);
760  SFREE(amta);
761  SFREE(namta);
762  SFREE(iamforc);
763  SFREE(iamload);
764  SFREE(iamboun);
765  }
766 
767  if(ntrans > 0){
768  RENEW(trab,double,7*ntrans);
769  }
770  else{SFREE(trab);SFREE(inotr);}
771 
772  if(ithermal[0] == 0){SFREE(t0);SFREE(t1);SFREE(t0g);SFREE(t1g);}
773  if((ithermal[0] == 0)||(nam<=0)){SFREE(iamt1);}
774 
775  if(ncs_>0){
776  RENEW(ics,ITG,ncs_);
777  SFREE(dcs);
778  }else if(npt_>0){SFREE(ics);}
779 
780  if(mcs>0){
781  RENEW(cs,double,17*mcs);
782  }else{
783  SFREE(cs);
784  }
785 
786  }else{
787 
788  /* reallocating space in all but the first step (>1) */
789 
790  RENEW(vold,double,mt*nk);
791 
792  /* if the SPC boundary conditions were changed in the present step,
793  they have to be rematched with those in the last step. Removed SPC
794  boundary conditions do not appear any more (this is different from
795  forces and loads, where removed forces or loads are reset to zero;
796  a removed SPC constraint does not have a numerical value any more) */
797 
798  NNEW(reorder,double,nboun);
799  NNEW(nreorder,ITG,nboun);
800  if(nbounold<nboun){
801  RENEW(xbounold,double,nboun);
802  RENEW(nodebounold,ITG,nboun);
803  RENEW(ndirbounold,ITG,nboun);
804  }
805  FORTRAN(spcmatch,(xboun,nodeboun,ndirboun,&nboun,xbounold,nodebounold,
806  ndirbounold,&nbounold,ikboun,ilboun,vold,reorder,nreorder,
807  mi));
808  RENEW(xbounold,double,nboun);
809  RENEW(nodebounold,ITG,nboun);
810  RENEW(ndirbounold,ITG,nboun);
811  SFREE(reorder); SFREE(nreorder);
812 
813  /* for additional forces or loads in the present step, the
814  corresponding slots in the force and load fields of the
815  previous steps are initialized */
816 
817  RENEW(xforcold,double,nforc);
818  for(i=nforcold;i<nforc;i++) xforcold[i]=0;
819 
820  RENEW(xloadold,double,2*nload);
821  for(i=2*nloadold;i<2*nload;i++) xloadold[i]=0;
822 
823  if(ithermal[0]!=0){
824  RENEW(t1old,double,nk);
825  }
826 
827  if(nam > 0) {
828  RENEW(amname,char,80*nam);
829  RENEW(namta,ITG,3*nam);
830  RENEW(amta,double,2*namta[3*nam-2]);
831  }
832 
833  }
834 
835  /* reallocating fields for all steps (>=1) */
836 
837  RENEW(co,double,3*nk);
838 
839  RENEW(nodeboun,ITG,nboun);
840  RENEW(ndirboun,ITG,nboun);
841  RENEW(typeboun,char,nboun+1);
842  RENEW(xboun,double,nboun);
843  RENEW(ikboun,ITG,nboun);
844  RENEW(ilboun,ITG,nboun);
845 
846  RENEW(nodeforc,ITG,2*nforc);
847  RENEW(ndirforc,ITG,nforc);
848  RENEW(xforc,double,nforc);
849  RENEW(ikforc,ITG,nforc);
850  RENEW(ilforc,ITG,nforc);
851 
852  /* temperature loading */
853 
854  if(ithermal[0] != 0){
855  RENEW(t0,double,nk);
856  RENEW(t1,double,nk);
857  if((ne1d!=0)||(ne2d!=0)){
858  RENEW(t0g,double,2*nk);
859  RENEW(t1g,double,2*nk);
860  }
861  if(nam > 0) {RENEW(iamt1,ITG,nk);}
862  }
863 
864  RENEW(nelemload,ITG,2*nload);
865  RENEW(sideload,char,20*nload);
866  RENEW(xload,double,2*nload);
867 
868  RENEW(cbody,char,81*nbody);
869  RENEW(ibody,ITG,3*nbody);
870  RENEW(xbody,double,7*nbody);
871  RENEW(xbodyold,double,7*nbody);
872 
873  RENEW(ipompc,ITG,nmpc);
874  RENEW(labmpc,char,20*nmpc+1);
875  RENEW(ikmpc,ITG,nmpc);
876  RENEW(ilmpc,ITG,nmpc);
877  RENEW(fmpc,double,nmpc);
878 
879  /* energy */
880 
881  if((nener==1)&&(nenerold==0)){
882  NNEW(ener,double,mi[0]*ne*2);
883  if((istep>1)&&(iperturb[0]>1)){
884  printf("*ERROR in CalculiX: in nonlinear calculations");
885  printf(" energy output must be selected in the first step");
886  FORTRAN(stop,());
887  }
888  }
889 
890  /* initial velocities and accelerations */
891 
892  if((nmethod==4)||(nmethod==5)||(nmethod==8)||(nmethod==9)||
893  ((abs(nmethod)==1)&&(iperturb[0]>=2))){
894  RENEW(veold,double,mt*nk);
895  }
896  else {SFREE(veold);}
897 
898  if((nmethod == 4)&&(iperturb[0]>1)) {
899  NNEW(accold,double,mt*nk);
900  }
901 
902  if(nam > 0) {
903  RENEW(iamforc,ITG,nforc);
904  RENEW(iamload,ITG,2*nload);
905  RENEW(iamboun,ITG,nboun);
906  }
907 
908  /* generate force convection elements */
909 
910  if(network==1){
911  ne0=ne;nkon0=nkon;nload1=nload;
912  RENEW(ipkon,ITG,ne+nload);
913  RENEW(lakon,char,8*(ne+nload));
914  RENEW(kon,ITG,nkon+9*nload);
915  NNEW(inodesd,ITG,nk);
916  RENEW(nelemload,ITG,4*nload);
917  RENEW(sideload,char,40*nload);
918 
919  FORTRAN(genadvecelem,(inodesd,ipkon,&ne,lakon,kon,&nload,
920  sideload,nelemload,&nkon));
921 
922  SFREE(inodesd);
923  RENEW(ipkon,ITG,ne);
924  RENEW(lakon,char,8*ne);
925  RENEW(kon,ITG,nkon);
926  RENEW(sti,double,6*mi[0]*ne);
927  RENEW(eme,double,6*mi[0]*ne);
928  if(iprestr>0) RENEW(prestr,double,6*mi[0]*ne);
929  if(nprop>0) RENEW(ielprop,ITG,ne);
930  if((ne1d!=0)||(ne2d!=0)) RENEW(offset,double,2*ne);
931  RENEW(nelemload,ITG,2*nload);
932  RENEW(sideload,char,20*nload);
933  RENEW(xload,double,2*nload);
934  RENEW(xloadold,double,2*nload);
935  if(nam>0){
936  RENEW(iamload,ITG,2*nload);
937  for(i=2*nload1;i<2*nload;i++)iamload[i]=0;
938  }
939  if(nener==1)RENEW(ener,double,mi[0]*ne*2);
940  if(norien>0)RENEW(ielorien,ITG,mi[2]*ne);
941  RENEW(ielmat,ITG,mi[2]*ne);
942  for(i=mi[2]*ne0;i<mi[2]*ne;i++)ielmat[i]=1;
943  }
944 
945  if(ntrans > 0){
946  RENEW(inotr,ITG,2*nk);
947  }
948 
949  /* calling the user routine ufaceload (can be empty) */
950 
951  if(ithermal[1]>=2){
952  NNEW(sideloadtemp,char,20*nload);
953  for(i=0;i<nload;i++){
954  strcpy1(&sideloadtemp[20*i],&sideload[20*i],20);
955  if((strcmp1(&sideload[20*i]," ")==0)&&
956  (strcmp1(&sideload[20*i+1]," ")!=0)){
957  strcpy1(&sideloadtemp[20*i],"F",1);
958  }
959  }
960  FORTRAN(ufaceload,(co,ipkon,kon,lakon,&nboun,nodeboun,
961  nelemload,sideloadtemp,&nload,&ne,&nk));
962  SFREE(sideloadtemp);
963  }
964 
965  /* decascading MPC's only necessary if MPC's changed */
966 
967  if(((istep == 1)||(ntrans>0)||(mpcend<0)||(nk!=nkold)||(nmpc!=nmpcold))&&(icascade==0)) {
968 
969  /* decascading the MPC's */
970 
971  printf(" Decascading the MPC's\n\n");
972 
973  callfrommain=1;
974  cascade(ipompc,&coefmpc,&nodempc,&nmpc,
975  &mpcfree,nodeboun,ndirboun,&nboun,ikmpc,
976  ilmpc,ikboun,ilboun,&mpcend,&mpcmult,
977  labmpc,&nk,&memmpc_,&icascade,&maxlenmpc,
978  &callfrommain,iperturb,ithermal);
979  }
980 
981  /* determining the matrix structure: changes if SPC's have changed */
982 
983  if((icascade==0)&&(nmethod<8)) printf(" Determining the structure of the matrix:\n");
984 
985  NNEW(nactdof,ITG,mt*nk);
986  NNEW(mast1,ITG,nzs[1]);
987  NNEW(irow,ITG,nzs[1]);
988 
989  if((mcs==0)||(cs[1]<0)){
990 
991  NNEW(icol,ITG,mt*nk);
992  NNEW(jq,ITG,mt*nk+1);
993  NNEW(ipointer,ITG,mt*nk);
994 
995  if((icascade==0)&&((nmethod<8)||(nmethod==11))){
996  if(nmethod==11){nmethodl=2;}else{nmethodl=nmethod;}
997  mastruct(&nk,kon,ipkon,lakon,&ne,nodeboun,ndirboun,&nboun,ipompc,
998  nodempc,&nmpc,nactdof,icol,jq,&mast1,&irow,&isolver,neq,
999  ikmpc,ilmpc,ipointer,nzs,&nmethodl,ithermal,
1000  ikboun,ilboun,iperturb,mi,&mortar,typeboun,labmpc);
1001  }
1002  else{neq[0]=1;neq[1]=1;neq[2]=1;}
1003  }
1004  else{
1005 
1006  NNEW(icol,ITG,8*nk);
1007  NNEW(jq,ITG,8*nk+1);
1008  NNEW(ipointer,ITG,8*nk);
1009 
1010  mastructcs(&nk,kon,ipkon,lakon,&ne,nodeboun,ndirboun,&nboun,
1011  ipompc,nodempc,&nmpc,nactdof,icol,jq,&mast1,&irow,&isolver,
1012  neq,ikmpc,ilmpc,ipointer,nzs,&nmethod,
1013  ics,cs,labmpc,&mcs,mi,&mortar);
1014  }
1015 
1016  SFREE(ipointer);SFREE(mast1);
1017  if((icascade==0)&&(nmethod<8))RENEW(irow,ITG,nzs[2]);
1018 
1019  /* nmethod=1: static analysis */
1020  /* nmethod=2: frequency analysis */
1021  /* nmethod=3: buckling analysis */
1022  /* nmethod=4: (linear or nonlinear) dynamic analysis */
1023  /* nmethod=5: steady state dynamics analysis */
1024  /* nmethod=6: Coriolis frequency calculation */
1025  /* nmethod=7: flutter frequency calculation */
1026 
1027 
1028  if((nmethod<=1)||(nmethod==11)||((iperturb[0]>1)&&(nmethod<8)))
1029  {
1030  if(iperturb[0]<2){
1031 
1032  linstatic(co,&nk,kon,ipkon,lakon,&ne,nodeboun,ndirboun,xboun,&nboun,
1033  ipompc,nodempc,coefmpc,labmpc,&nmpc,nodeforc,ndirforc,xforc,
1034  &nforc, nelemload,sideload,xload,&nload,
1035  nactdof,&icol,jq,&irow,neq,&nzl,&nmethod,ikmpc,
1036  ilmpc,ikboun,ilboun,elcon,nelcon,rhcon,nrhcon,
1037  alcon,nalcon,alzero,ielmat,ielorien,&norien,orab,&ntmat_,
1038  t0,t1,t1old,ithermal,prestr,&iprestr, vold,iperturb,sti,nzs,
1039  &kode,filab,eme,&iexpl,plicon,
1040  nplicon,plkcon,nplkcon,xstate,&npmat_,matname,
1041  &isolver,mi,&ncmat_,&nstate_,cs,&mcs,&nkon,ener,
1042  xbounold,xforcold,xloadold,amname,amta,namta,
1043  &nam,iamforc,iamload,iamt1,iamboun,&ttime,
1044  output,set,&nset,istartset,iendset,ialset,&nprint,prlab,
1045  prset,&nener,trab,inotr,&ntrans,fmpc,cbody,ibody,xbody,&nbody,
1046  xbodyold,timepar,thicke,jobnamec,tieset,&ntie,&istep,&nmat,
1047  ielprop,prop,typeboun);
1048 
1049  }
1050 
1051  else{
1052 
1053  mpcinfo[0]=memmpc_;mpcinfo[1]=mpcfree;mpcinfo[2]=icascade;
1054  mpcinfo[3]=maxlenmpc;
1055 
1056  nonlingeo(&co,&nk,&kon,&ipkon,&lakon,&ne,nodeboun,ndirboun,xboun,&nboun,
1057  &ipompc,&nodempc,&coefmpc,&labmpc,&nmpc,nodeforc,ndirforc,xforc,
1058  &nforc, nelemload,sideload,xload,&nload,
1059  nactdof,&icol,jq,&irow,neq,&nzl,&nmethod,&ikmpc,
1060  &ilmpc,ikboun,ilboun,elcon,nelcon,rhcon,nrhcon,
1061  alcon,nalcon,alzero,&ielmat,&ielorien,&norien,orab,&ntmat_,
1062  t0,t1,t1old,ithermal,prestr,&iprestr,
1063  &vold,iperturb,sti,nzs,&kode,filab,&idrct,jmax,
1064  jout,timepar,eme,xbounold,xforcold,xloadold,
1065  veold,accold,amname,amta,namta,
1066  &nam,iamforc,iamload,iamt1,&alpha,
1067  &iexpl,iamboun,plicon,nplicon,plkcon,nplkcon,
1068  &xstate,&npmat_,&istep,&ttime,matname,qaold,mi,
1069  &isolver,&ncmat_,&nstate_,&iumat,cs,&mcs,&nkon,&ener,
1070  mpcinfo,output,
1071  shcon,nshcon,cocon,ncocon,physcon,&nflow,ctrl,
1072  set,&nset,istartset,iendset,ialset,&nprint,prlab,
1073  prset,&nener,ikforc,ilforc,trab,inotr,&ntrans,&fmpc,
1074  cbody,ibody,xbody,&nbody,xbodyold,ielprop,prop,
1075  &ntie,tieset,&itpamp,&iviewfile,jobnamec,tietol,&nslavs,thicke,
1076  ics,&nintpoint,&mortar,
1077  &ifacecount,typeboun,&islavsurf,&pslavsurf,&clearini,&nmat,
1078  xmodal,&iaxial,&inext);
1079 
1080  memmpc_=mpcinfo[0];mpcfree=mpcinfo[1];icascade=mpcinfo[2];
1081  maxlenmpc=mpcinfo[3];
1082 
1083  }
1084  }else if(nmethod==2){
1085 
1086  /* FREQUENCY ANALYSIS */
1087 
1088  if((mcs==0)||(cs[1]<0)){
1089 #ifdef ARPACK
1090 
1091  mpcinfo[0]=memmpc_;mpcinfo[1]=mpcfree;mpcinfo[2]=icascade;
1092  mpcinfo[3]=maxlenmpc;
1093 
1094  arpack(co,&nk,&kon,&ipkon,&lakon,&ne,nodeboun,ndirboun,xboun,&nboun,
1095  ipompc,nodempc,coefmpc,labmpc,&nmpc,nodeforc,ndirforc,xforc,
1096  &nforc, nelemload,sideload,xload,&nload,
1097  nactdof,icol,jq,&irow,neq,&nzl,&nmethod,ikmpc,
1098  ilmpc,ikboun,ilboun,elcon,nelcon,rhcon,nrhcon,
1099  shcon,nshcon,cocon,ncocon,
1100  alcon,nalcon,alzero,&ielmat,&ielorien,&norien,orab,&ntmat_,
1101  t0,t1,t1old,ithermal,prestr,&iprestr,vold,iperturb,sti,nzs,
1102  &kode,mei,fei,filab,
1103  eme,&iexpl,plicon,nplicon,plkcon,nplkcon,
1104  &xstate,&npmat_,matname,mi,&ncmat_,&nstate_,&ener,jobnamec,
1105  output,set,&nset,istartset,iendset,ialset,&nprint,prlab,
1106  prset,&nener,&isolver,trab,inotr,&ntrans,&ttime,fmpc,cbody,
1107  ibody,xbody,&nbody,thicke,&nslavs,tietol,&nkon,mpcinfo,
1108  &ntie,&istep,&mcs,ics,tieset,cs,&nintpoint,&mortar,&ifacecount,
1109  &islavsurf,&pslavsurf,&clearini,&nmat,typeboun,ielprop,prop);
1110 
1111  memmpc_=mpcinfo[0];mpcfree=mpcinfo[1];icascade=mpcinfo[2];
1112  maxlenmpc=mpcinfo[3];
1113 
1114 #else
1115  printf("*ERROR in CalculiX: the ARPACK library is not linked\n\n");
1116  FORTRAN(stop,());
1117 #endif
1118 
1119  }else{
1120 
1121 #ifdef ARPACK
1122 
1123  mpcinfo[0]=memmpc_;mpcinfo[1]=mpcfree;mpcinfo[2]=icascade;
1124  mpcinfo[3]=maxlenmpc;
1125 
1126  arpackcs(co,&nk,&kon,&ipkon,&lakon,&ne,nodeboun,ndirboun,
1127  xboun,&nboun,
1128  ipompc,nodempc,coefmpc,labmpc,&nmpc,nodeforc,ndirforc,xforc,
1129  &nforc, nelemload,sideload,xload,&nload,
1130  nactdof,icol,jq,&irow,neq,&nzl,&nmethod,ikmpc,
1131  ilmpc,ikboun,ilboun,elcon,nelcon,rhcon,nrhcon,
1132  alcon,nalcon,alzero,&ielmat,&ielorien,&norien,orab,&ntmat_,
1133  t0,t1,t1old,ithermal,prestr,&iprestr,
1134  vold,iperturb,sti,nzs,&kode,mei,fei,filab,
1135  eme,&iexpl,plicon,nplicon,plkcon,nplkcon,
1136  &xstate,&npmat_,matname,mi,ics,cs,&mpcend,&ncmat_,
1137  &nstate_,&mcs,&nkon,&ener,jobnamec,output,set,&nset,istartset,
1138  iendset,ialset,&nprint,prlab,
1139  prset,&nener,&isolver,trab,inotr,&ntrans,&ttime,fmpc,cbody,
1140  ibody,xbody,&nbody,&nevtot,thicke,&nslavs,tietol,mpcinfo,
1141  &ntie,&istep,tieset,&nintpoint,&mortar,&ifacecount,&islavsurf,
1142  &pslavsurf,&clearini,&nmat,typeboun,ielprop,prop);
1143 
1144  memmpc_=mpcinfo[0];mpcfree=mpcinfo[1];icascade=mpcinfo[2];
1145  maxlenmpc=mpcinfo[3];
1146 
1147 #else
1148  printf("*ERROR in CalculiX: the ARPACK library is not linked\n\n");
1149  FORTRAN(stop,());
1150 #endif
1151 
1152  }
1153  }else if(nmethod==3){
1154 
1155 #ifdef ARPACK
1156  arpackbu(co,&nk,kon,ipkon,lakon,&ne,nodeboun,ndirboun,xboun,&nboun,
1157  ipompc,nodempc,coefmpc,labmpc,&nmpc,nodeforc,ndirforc,xforc,
1158  &nforc,
1159  nelemload,sideload,xload,&nload,
1160  nactdof,icol,jq,irow,neq,&nzl,&nmethod,ikmpc,
1161  ilmpc,ikboun,ilboun,elcon,nelcon,rhcon,nrhcon,
1162  alcon,nalcon,alzero,ielmat,ielorien,&norien,orab,&ntmat_,
1163  t0,t1,t1old,ithermal,prestr,&iprestr,
1164  vold,iperturb,sti,nzs,&kode,mei,fei,filab,
1165  eme,&iexpl,plicon,nplicon,plkcon,nplkcon,
1166  xstate,&npmat_,matname,mi,&ncmat_,&nstate_,ener,output,
1167  set,&nset,istartset,iendset,ialset,&nprint,prlab,
1168  prset,&nener,&isolver,trab,inotr,&ntrans,&ttime,fmpc,cbody,
1169  ibody,xbody,&nbody,thicke,jobnamec,&nmat,ielprop,prop);
1170 #else
1171  printf("*ERROR in CalculiX: the ARPACK library is not linked\n\n");
1172  FORTRAN(stop,());
1173 #endif
1174  }
1175  else if(nmethod==4)
1176  {
1177  if((ne1d!=0)||(ne2d!=0)){
1178  printf(" *WARNING: 1-D or 2-D elements may cause problems in modal dynamic calculations\n");
1179  printf(" ensure that point loads defined in a *MODAL DYNAMIC step\n");
1180  printf(" and applied to nodes belonging to 1-D or 2-D elements have been\n");
1181  printf(" applied to the same nodes in the preceding FREQUENCY step with\n");
1182  printf(" magnitude zero; look at example shellf.inp for a guideline.\n\n");}
1183 
1184  printf(" Composing the dynamic response from the eigenmodes\n\n");
1185 
1186  dyna(&co,&nk,&kon,&ipkon,&lakon,&ne,&nodeboun,&ndirboun,&xboun,&nboun,
1187  &ipompc,&nodempc,&coefmpc,&labmpc,&nmpc,nodeforc,ndirforc,xforc,&nforc,
1188  nelemload,sideload,xload,&nload,
1189  &nactdof,neq,&nzl,icol,irow,&nmethod,&ikmpc,&ilmpc,&ikboun,&ilboun,
1190  elcon,nelcon,rhcon,nrhcon,cocon,ncocon,
1191  alcon,nalcon,alzero,&ielmat,&ielorien,&norien,orab,&ntmat_,&t0,
1192  &t1,ithermal,prestr,&iprestr,&vold,iperturb,&sti,nzs,
1193  timepar,xmodal,&veold,amname,amta,
1194  namta,&nam,iamforc,iamload,&iamt1,
1195  jout,&kode,filab,&eme,xforcold,xloadold,
1196  &t1old,&iamboun,&xbounold,&iexpl,plicon,
1197  nplicon,plkcon,nplkcon,&xstate,&npmat_,matname,
1198  mi,&ncmat_,&nstate_,&ener,jobnamec,&ttime,set,&nset,
1199  istartset,iendset,&ialset,&nprint,prlab,
1200  prset,&nener,trab,&inotr,&ntrans,&fmpc,cbody,ibody,xbody,&nbody,
1201  xbodyold,&istep,&isolver,jq,output,&mcs,&nkon,&mpcend,ics,cs,
1202  &ntie,tieset,&idrct,jmax,ctrl,&itpamp,tietol,&nalset,
1203  ikforc,ilforc,thicke,&nslavs,&nmat,typeboun,ielprop,prop);
1204  }
1205  else if(nmethod==5)
1206  {
1207  if((ne1d!=0)||(ne2d!=0)){
1208  printf(" *WARNING: 1-D or 2-D elements may cause problems in steady state calculations\n");
1209  printf(" ensure that point loads defined in a *STEADY STATE DYNAMICS step\n");
1210  printf(" and applied to nodes belonging to 1-D or 2-D elements have been\n");
1211  printf(" applied to the same nodes in the preceding FREQUENCY step with\n");
1212  printf(" magnitude zero; look at example shellf.inp for a guideline.\n\n");}
1213 
1214  printf(" Composing the steady state response from the eigenmodes\n\n");
1215 
1216  steadystate(&co,&nk,&kon,&ipkon,&lakon,&ne,&nodeboun,&ndirboun,&xboun,&nboun,
1217  &ipompc,&nodempc,&coefmpc,&labmpc,&nmpc,nodeforc,ndirforc,xforc,&nforc,
1218  nelemload,sideload,xload,&nload,
1219  &nactdof,neq,&nzl,icol,irow,&nmethod,&ikmpc,&ilmpc,&ikboun,&ilboun,
1220  elcon,nelcon,rhcon,nrhcon,cocon,ncocon,
1221  alcon,nalcon,alzero,&ielmat,&ielorien,&norien,orab,&ntmat_,&t0,
1222  &t1,ithermal,prestr,&iprestr,&vold,iperturb,sti,nzs,
1223  timepar,xmodal,&veold,amname,amta,
1224  namta,&nam,iamforc,iamload,&iamt1,
1225  jout,&kode,filab,&eme,xforcold,xloadold,
1226  &t1old,&iamboun,&xbounold,&iexpl,plicon,
1227  nplicon,plkcon,nplkcon,xstate,&npmat_,matname,
1228  mi,&ncmat_,&nstate_,&ener,jobnamec,&ttime,set,&nset,
1229  istartset,iendset,ialset,&nprint,prlab,
1230  prset,&nener,trab,&inotr,&ntrans,&fmpc,cbody,ibody,xbody,&nbody,
1231  xbodyold,&istep,&isolver,jq,output,&mcs,&nkon,ics,cs,&mpcend,
1232  ctrl,ikforc,ilforc,thicke,&nmat,typeboun,ielprop,prop);
1233  }
1234  else if((nmethod==6)||(nmethod==7))
1235  {
1236 
1237  printf(" Composing the complex eigenmodes from the real eigenmodes\n\n");
1238 
1239  complexfreq(&co,&nk,&kon,&ipkon,&lakon,&ne,&nodeboun,&ndirboun,&xboun,&nboun,
1240  &ipompc,&nodempc,&coefmpc,&labmpc,&nmpc,nodeforc,ndirforc,xforc,&nforc,
1241  nelemload,sideload,xload,&nload,
1242  &nactdof,neq,&nzl,icol,irow,&nmethod,&ikmpc,&ilmpc,&ikboun,&ilboun,
1243  elcon,nelcon,rhcon,nrhcon,cocon,ncocon,
1244  alcon,nalcon,alzero,&ielmat,&ielorien,&norien,orab,&ntmat_,&t0,
1245  &t1,ithermal,prestr,&iprestr,&vold,iperturb,&sti,nzs,
1246  timepar,xmodal,&veold,amname,amta,
1247  namta,&nam,iamforc,iamload,&iamt1,
1248  jout,&kode,filab,&eme,xforcold,xloadold,
1249  &t1old,&iamboun,&xbounold,&iexpl,plicon,
1250  nplicon,plkcon,nplkcon,xstate,&npmat_,matname,
1251  mi,&ncmat_,&nstate_,&ener,jobnamec,&ttime,set,&nset,
1252  istartset,iendset,&ialset,&nprint,prlab,
1253  prset,&nener,trab,&inotr,&ntrans,&fmpc,cbody,ibody,xbody,&nbody,
1254  xbodyold,&istep,&isolver,jq,output,&mcs,&nkon,&mpcend,ics,cs,
1255  &ntie,tieset,&idrct,jmax,ctrl,&itpamp,tietol,&nalset,
1256  ikforc,ilforc,thicke,jobnamef,mei,&nmat,ielprop,prop);
1257  }
1258  else if(nmethod>7){
1259 
1260  mpcinfo[0]=memmpc_;mpcinfo[1]=mpcfree;mpcinfo[2]=icascade;
1261  mpcinfo[3]=maxlenmpc;
1262 
1263  electromagnetics(&co,&nk,&kon,&ipkon,&lakon,&ne,nodeboun,
1264  ndirboun,xboun,&nboun,
1265  &ipompc,&nodempc,&coefmpc,&labmpc,&nmpc,nodeforc,ndirforc,xforc,
1266  &nforc,&nelemload,&sideload,xload,&nload,
1267  nactdof,&icol,jq,&irow,neq,&nzl,&nmethod,&ikmpc,
1268  &ilmpc,ikboun,ilboun,elcon,nelcon,rhcon,nrhcon,
1269  alcon,nalcon,alzero,&ielmat,&ielorien,&norien,orab,&ntmat_,
1270  t0,t1,t1old,ithermal,prestr,&iprestr,
1271  &vold,iperturb,sti,nzs,&kode,filab,&idrct,jmax,
1272  jout,timepar,eme,xbounold,xforcold,xloadold,
1273  veold,accold,amname,amta,namta,
1274  &nam,iamforc,&iamload,iamt1,&alpha,
1275  &iexpl,iamboun,plicon,nplicon,plkcon,nplkcon,
1276  &xstate,&npmat_,&istep,&ttime,matname,qaold,mi,
1277  &isolver,&ncmat_,&nstate_,&iumat,cs,&mcs,&nkon,&ener,
1278  mpcinfo,output,
1279  shcon,nshcon,cocon,ncocon,physcon,&nflow,ctrl,
1280  &set,&nset,&istartset,&iendset,&ialset,&nprint,prlab,
1281  prset,&nener,ikforc,ilforc,trab,inotr,&ntrans,&fmpc,
1282  cbody,ibody,xbody,&nbody,xbodyold,ielprop,prop,
1283  &ntie,&tieset,&itpamp,&iviewfile,jobnamec,&tietol,&nslavs,thicke,
1284  ics,&nalset,&nmpc_,&nmat,typeboun,&iaxial,&idefload,&nload_);
1285 
1286  memmpc_=mpcinfo[0];mpcfree=mpcinfo[1];icascade=mpcinfo[2];
1287  maxlenmpc=mpcinfo[3];
1288  }
1289 
1290  SFREE(nactdof);
1291  SFREE(icol);
1292  SFREE(jq);
1293  SFREE(irow);
1294 
1295  /* deleting the perturbation loads and temperatures */
1296 
1297  if((iperturb[0] == 1)&&(nmethod==3)) {
1298  nforc=0;
1299  nload=0;
1300  nbody=0;
1301  if(ithermal[0] == 1) {
1302  for(k=0;k<nk;++k){
1303  t1[k]=t0[k];
1304  }
1305  }
1306  }else{
1307  nbounold=nboun;
1308  for (i=0;i<nboun;i++) {
1309  nodebounold[i]=nodeboun[i];
1310  ndirbounold[i]=ndirboun[i];
1311  }
1312  nforcold=nforc;
1313  nloadold=nload;
1314  nbodyold=nbody;
1315 
1316  /* resetting the amplitude to none except for time=total time amplitudes */
1317 
1318  if(nam > 0) {
1319  for (i=0;i<nboun;i++) {
1320  if(iamboun[i]>0){
1321  if(namta[3*iamboun[i]-1]>0){
1322  iamboun[i]=0;
1323  xboun[i]=xbounold[i];}
1324  }
1325  }
1326  for (i=0;i<nforc;i++){
1327  if(iamforc[i]>0){
1328  if(namta[3*iamforc[i]-1]>0){
1329  iamforc[i]=0;
1330  xforc[i]=xforcold[i];}
1331  }
1332  }
1333  for (i=0;i<2*nload;i++){
1334  if(iamload[i]>0){
1335  if(namta[3*iamload[i]-1]>0){
1336  iamload[i]=0;
1337  xload[i]=xloadold[i];}
1338  }
1339  }
1340  for (i=1;i<3*nbody;i=i+3){
1341  if(ibody[i]>0){
1342  if(namta[3*ibody[i]-1]>0){
1343  ibody[i]=0;
1344  xbody[7*(i-1)/3]=xbodyold[7*(i-1)/3];}
1345  }
1346  }
1347  if(ithermal[0]==1) {
1348  if(iamt1[i]>0){
1349  if(namta[3*iamt1[i]-1]>0){
1350  iamt1[i]=0;
1351  t1[i]=t1old[i];}
1352  }
1353  }
1354  }
1355  }
1356 
1357  /* removing the advective elements, if any */
1358 
1359  if(network==1){
1360  ne=ne0;nkon=nkon0;
1361  RENEW(ipkon,ITG,ne);
1362  RENEW(lakon,char,8*ne);
1363  RENEW(kon,ITG,nkon);
1364  RENEW(sti,double,6*mi[0]*ne);
1365  RENEW(eme,double,6*mi[0]*ne);
1366  if(iprestr>0) RENEW(prestr,double,6*mi[0]*ne);
1367  if(nprop>0) RENEW(ielprop,ITG,ne);
1368  if((ne1d!=0)||(ne2d!=0)) RENEW(offset,double,2*ne);
1369  if(nener==1)RENEW(ener,double,mi[0]*ne*2);
1370  if(norien>0)RENEW(ielorien,ITG,mi[2]*ne);
1371  RENEW(ielmat,ITG,mi[2]*ne);
1372 
1373  /* reactivating the original load labels */
1374 
1375  for(i=nload-1;i>=nload0;i--){
1376  if(strcmp2(&sideload[20*i]," ",20)==0){
1377  iload=nelemload[2*i+1];
1378  strcpy1(&sideload[20*(iload-1)],"F",1);
1379  }
1380  }
1381 
1382  }
1383 
1384  nload=nload0;
1385 
1386  if((nmethod == 4)&&(iperturb[0]>1)) SFREE(accold);
1387 
1388  if(irstrt>0){
1389  jrstrt++;
1390  if(jrstrt==irstrt){
1391  jrstrt=0;
1392  FORTRAN(restartwrite,(&istep,&nset,&nload,&nforc,&nboun,&nk,&ne,
1393  &nmpc,&nalset,&nmat,&ntmat_,&npmat_,&norien,&nam,&nprint,
1394  mi,&ntrans,&ncs_,&namtot_,&ncmat_,&mpcend,&maxlenmpc,&ne1d,
1395  &ne2d,&nflow,&nlabel,&iplas,&nkon,ithermal,&nmethod,iperturb,
1396  &nstate_,&nener,set,istartset,iendset,ialset,co,kon,ipkon,
1397  lakon,nodeboun,ndirboun,iamboun,xboun,ikboun,ilboun,ipompc,
1398  nodempc,coefmpc,labmpc,ikmpc,ilmpc,nodeforc,ndirforc,iamforc,
1399  xforc,ikforc,ilforc,nelemload,iamload,sideload,xload,
1400  elcon,nelcon,rhcon,nrhcon,alcon,nalcon,
1401  alzero,plicon,nplicon,plkcon,nplkcon,orname,orab,ielorien,
1402  trab,inotr,amname,amta,namta,t0,t1,iamt1,veold,
1403  ielmat,matname,prlab,prset,filab,vold,nodebounold,
1404  ndirbounold,xbounold,xforcold,xloadold,t1old,eme,
1405  iponor,xnor,knor,thicke,offset,iponoel,inoel,rig,
1406  shcon,nshcon,cocon,ncocon,ics,
1407  sti,ener,xstate,jobnamec,infree,prestr,&iprestr,cbody,
1408  ibody,xbody,&nbody,xbodyold,&ttime,qaold,cs,&mcs,output,
1409  physcon,ctrl,typeboun,fmpc,tieset,&ntie,tietol,&nslavs,t0g,t1g,
1410  &nprop,ielprop,prop,&mortar,&nintpoint,&ifacecount,islavsurf,
1411  pslavsurf,clearini));
1412  }
1413  }
1414 
1415 }
1416 
1417 FORTRAN(closefile,());
1418 
1419 strcpy(fneig,jobnamec);
1420 strcat(fneig,".frd");
1421 if((f1=fopen(fneig,"ab"))==NULL){
1422  printf("*ERROR in frd: cannot open frd file for writing...");
1423  exit(0);
1424 }
1425 #ifdef EXODUSII
1426  if(strcmp1(output,"exo")==1){
1427 #endif
1428 fprintf(f1," 9999\n");
1429 fclose(f1);
1430 #ifdef EXODUSII
1431  }
1432 #endif
1433 
1434 /* deallocating the fields
1435  this section is addressed immediately after leaving calinput */
1436 
1437 SFREE(ipoinpc);SFREE(inpc);SFREE(inp);
1438 
1439 if(ncs_>0) SFREE(ics);
1440 if((ncs_<=0)&&(npt_>0)) SFREE(dcs);
1441 if(mcs>0) SFREE(cs);
1442 SFREE(tieset);SFREE(tietol);
1443 
1444 SFREE(co);SFREE(kon);SFREE(ipkon);SFREE(lakon);
1445 
1446 SFREE(nodeboun);SFREE(ndirboun);SFREE(typeboun);SFREE(xboun);SFREE(ikboun);
1447 SFREE(ilboun);SFREE(nodebounold);SFREE(ndirbounold);SFREE(xbounold);
1448 
1449 SFREE(ipompc);SFREE(labmpc);SFREE(ikmpc);SFREE(ilmpc);SFREE(fmpc);
1450 SFREE(nodempc);SFREE(coefmpc);
1451 
1452 SFREE(nodeforc);SFREE(ndirforc);SFREE(xforc);SFREE(ikforc);SFREE(ilforc);
1453 SFREE(xforcold);
1454 
1455 SFREE(nelemload);SFREE(sideload);SFREE(xload);SFREE(xloadold);
1456 
1457 SFREE(cbody);SFREE(ibody);SFREE(xbody);SFREE(xbodyold);
1458 
1459 if(nam>0){SFREE(iamboun);SFREE(iamforc);SFREE(iamload);SFREE(amname);
1460  SFREE(amta);SFREE(namta);}
1461 
1462 SFREE(set);SFREE(istartset);SFREE(iendset);SFREE(ialset);
1463 
1464 SFREE(elcon);SFREE(nelcon);SFREE(rhcon);SFREE(nrhcon);SFREE(shcon);SFREE(nshcon);
1465 SFREE(cocon);SFREE(ncocon);SFREE(alcon);SFREE(nalcon);SFREE(alzero);
1466 if(nprop>0){SFREE(ielprop);SFREE(prop);}
1467 if(npmat_>0){SFREE(plicon);SFREE(nplicon);SFREE(plkcon);SFREE(nplkcon);}
1468 
1469 if(norien>0){SFREE(orname);SFREE(orab);SFREE(ielorien);}
1470 if(ntrans>0){SFREE(trab);SFREE(inotr);}
1471 if(iprestr>0){SFREE(prestr);}
1472 
1473 if(ithermal[0]!=0){
1474  SFREE(t0);SFREE(t1);SFREE(t1old);
1475  if(nam>0) SFREE(iamt1);
1476  if((ne1d!=0)||(ne2d!=0)){SFREE(t0g);SFREE(t1g);}
1477 }
1478 
1479 SFREE(prlab);SFREE(prset);SFREE(filab);SFREE(xmodal);
1480 
1481 SFREE(ielmat);SFREE(matname);
1482 
1483 SFREE(sti);SFREE(eme);SFREE(ener);SFREE(xstate);
1484 
1485 SFREE(vold);SFREE(veold);
1486 
1487 if((ne1d!=0)||(ne2d!=0)){
1488  SFREE(iponor);SFREE(xnor);SFREE(knor);SFREE(thicke);SFREE(offset);
1489  SFREE(iponoel);SFREE(inoel);SFREE(rig);
1490 }
1491 
1492 SFREE(islavsurf);
1493 if(mortar==1){SFREE(pslavsurf);SFREE(clearini);}
1494 
1495 #ifdef CALCULIX_MPI
1496 MPI_Finalize();
1497 #endif
1498 
1499  return 0;
1500 
1501 }
void cascade(ITG *ipompc, double **coefmpcp, ITG **nodempcp, ITG *nmpc, ITG *mpcfree, ITG *nodeboun, ITG *ndirboun, ITG *nboun, ITG *ikmpc, ITG *ilmpc, ITG *ikboun, ITG *ilboun, ITG *mpcend, ITG *mpcmult, char *labmpc, ITG *nk, ITG *memmpc_, ITG *icascade, ITG *maxlenmpc, ITG *callfrommain, ITG *iperturb, ITG *ithermal)
Definition: cascade.c:34
ITG strcmp2(const char *s1, const char *s2, ITG length)
Definition: strcmp2.c:24
subroutine ufaceload(co, ipkon, kon, lakon, nboun, nodeboun, nelemload, sideload, nload, ne, nk)
Definition: ufaceload.f:19
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))
void mastruct(ITG *nk, ITG *kon, ITG *ipkon, char *lakon, ITG *ne, ITG *nodeboun, ITG *ndirboun, ITG *nboun, ITG *ipompc, ITG *nodempc, ITG *nmpc, ITG *nactdof, ITG *icol, ITG *jq, ITG **mast1p, ITG **irowp, ITG *isolver, ITG *neq, ITG *ikmpc, ITG *ilmpc, ITG *ipointer, ITG *nzs, ITG *nmethod, ITG *ithermal, ITG *ikboun, ITG *ilboun, ITG *iperturb, ITG *mi, ITG *mortar, char *typeboun, char *labmpc)
Definition: mastruct.c:27
ITG strcpy1(char *s1, const char *s2, ITG length)
Definition: strcpy1.c:24
ITG strcmp1(const char *s1, const char *s2)
Definition: strcmp1.c:24
void readinput(char *jobnamec, char **inpcp, ITG *nline, ITG *nset, ITG *ipoinp, ITG **inpp, ITG **ipoinpcp, ITG *ithermal)
Definition: readinput.c:25
void arpackbu(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 *icol, ITG *jq, ITG *irow, 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, ITG *mei, double *fei, char *filab, double *eme, ITG *iexpl, double *plicon, ITG *nplicon, double *plkcon, ITG *nplkcon, double *xstate, ITG *npmat_, char *matname, ITG *mi, ITG *ncmat_, ITG *nstate_, double *ener, char *output, char *set, ITG *nset, ITG *istartset, ITG *iendset, ITG *ialset, ITG *nprint, char *prlab, char *prset, ITG *nener, ITG *isolver, double *trab, ITG *inotr, ITG *ntrans, double *ttime, double *fmpc, char *cbody, ITG *ibody, double *xbody, ITG *nbody, double *thicke, char *jobnamec, ITG *nmat, ITG *ielprop, double *prop)
subroutine openfile(jobname, output)
Definition: openfile.f:19
void tiedcontact(ITG *ntie, char *tieset, ITG *nset, char *set, ITG *istartset, ITG *iendset, ITG *ialset, char *lakon, ITG *ipkon, ITG *kon, double *tietol, ITG *nmpc, ITG *mpcfree, ITG *memmpc_, ITG **ipompcp, char **labmpcp, ITG **ikmpcp, ITG **ilmpcp, double **fmpcp, ITG **nodempcp, double **coefmpcp, ITG *ithermal, double *co, double *vold, ITG *cfd, ITG *nmpc_, ITG *mi, ITG *nk, ITG *istep, ITG *ikboun, ITG *nboun, char *kind1, char *kind2)
Definition: tiedcontact.c:23
#define DMEMSET(a, b, c, d)
Definition: CalculiX.h:45
void electromagnetics(double **co, ITG *nk, ITG **konp, ITG **ipkonp, char **lakonp, ITG *ne, ITG *nodeboun, ITG *ndirboun, double *xboun, ITG *nboun, ITG **ipompcp, ITG **nodempcp, double **coefmpcp, char **labmpcp, ITG *nmpc, ITG *nodeforc, ITG *ndirforc, double *xforc, ITG *nforc, ITG **nelemloadp, char **sideloadp, double *xload, ITG *nload, ITG *nactdof, ITG **icolp, ITG *jq, ITG **irowp, ITG *neq, ITG *nzl, ITG *nmethod, ITG **ikmpcp, ITG **ilmpcp, ITG *ikboun, ITG *ilboun, double *elcon, ITG *nelcon, double *rhcon, ITG *nrhcon, double *alcon, ITG *nalcon, double *alzero, ITG **ielmatp, ITG **ielorienp, 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, ITG *idrct, ITG *jmax, ITG *jout, double *timepar, double *eme, double *xbounold, double *xforcold, double *xloadold, double *veold, double *accold, char *amname, double *amta, ITG *namta, ITG *nam, ITG *iamforc, ITG **iamloadp, ITG *iamt1, double *alpha, ITG *iexpl, ITG *iamboun, double *plicon, ITG *nplicon, double *plkcon, ITG *nplkcon, double **xstatep, ITG *npmat_, ITG *istep, double *ttime, char *matname, double *qaold, ITG *mi, ITG *isolver, ITG *ncmat_, ITG *nstate_, ITG *iumat, double *cs, ITG *mcs, ITG *nkon, double **ener, ITG *mpcinfo, char *output, double *shcon, ITG *nshcon, double *cocon, ITG *ncocon, double *physcon, ITG *nflow, double *ctrl, char **setp, ITG *nset, ITG **istartsetp, ITG **iendsetp, ITG **ialsetp, ITG *nprint, char *prlab, char *prset, ITG *nener, ITG *ikforc, ITG *ilforc, double *trab, ITG *inotr, ITG *ntrans, double **fmpcp, char *cbody, ITG *ibody, double *xbody, ITG *nbody, double *xbodyold, ITG *ielprop, double *prop, ITG *ntie, char **tiesetp, ITG *itpamp, ITG *iviewfile, char *jobnamec, double **tietolp, ITG *nslavs, double *thicke, ITG *ics, ITG *nalset, ITG *nmpc_, ITG *nmat, char *typeboun, ITG *iaxial, ITG **idefloadp, ITG *nload_)
subroutine stop()
Definition: stop.f:19
subroutine restartwrite(istepnew, nset, nload, nforc, nboun, nk, ne, nmpc, nalset, nmat, ntmat_, npmat_, norien, nam, nprint, mi, ntrans, ncs_, namtot_, ncmat_, mpcend, maxlenmpc, ne1d, ne2d, nflow, nlabel, iplas, nkon, ithermal, nmethod, iperturb, nstate_, nener, set, istartset, iendset, ialset, co, kon, ipkon, lakon, nodeboun, ndirboun, iamboun, xboun, ikboun, ilboun, ipompc, nodempc, coefmpc, labmpc, ikmpc, ilmpc, nodeforc, ndirforc, iamforc, xforc, ikforc, ilforc, nelemload, iamload, sideload, xload, elcon, nelcon, rhcon, nrhcon, alcon, nalcon, alzero, plicon, nplicon, plkcon, nplkcon, orname, orab, ielorien, trab, inotr, amname, amta, namta, t0, t1, iamt1, veold, ielmat, matname, prlab, prset, filab, vold, nodebounold, ndirbounold, xbounold, xforcold, xloadold, t1old, eme, iponor, xnor, knor, thicke, offset, iponoel, inoel, rig, shcon, nshcon, cocon, ncocon, ics, sti, ener, xstate, jobnamec, infree, prestr, iprestr, cbody, ibody, xbody, nbody, xbodyold, ttime, qaold, cs, mcs, output, physcon, ctrl, typeboun, fmpc, tieset, ntie, tietol, nslavs, t0g, t1g, nprop, ielprop, prop, mortar, nintpoint, ifacecount, islavsurf, pslavsurf, clearini)
Definition: restartwrite.f:19
void dyna(double **cop, ITG *nk, ITG **konp, ITG **ipkonp, char **lakonp, ITG *ne, ITG **nodebounp, ITG **ndirbounp, double **xbounp, ITG *nboun, ITG **ipompcp, ITG **nodempcp, double **coefmpcp, char **labmpcp, ITG *nmpc, ITG *nodeforc, ITG *ndirforc, double *xforc, ITG *nforc, ITG *nelemload, char *sideload, double *xload, ITG *nload, ITG **nactdofp, ITG *neq, ITG *nzl, ITG *icol, ITG *irow, ITG *nmethod, ITG **ikmpcp, ITG **ilmpcp, ITG **ikbounp, ITG **ilbounp, double *elcon, ITG *nelcon, double *rhcon, ITG *nrhcon, double *cocon, ITG *ncocon, double *alcon, ITG *nalcon, double *alzero, ITG **ielmatp, ITG **ielorienp, ITG *norien, double *orab, ITG *ntmat_, double **t0p, double **t1p, ITG *ithermal, double *prestr, ITG *iprestr, double **voldp, ITG *iperturb, double **stip, ITG *nzs, double *timepar, double *xmodal, double **veoldp, char *amname, double *amta, ITG *namta, ITG *nam, ITG *iamforc, ITG *iamload, ITG **iamt1p, ITG *jout, ITG *kode, char *filab, double **emep, double *xforcold, double *xloadold, double **t1oldp, ITG **iambounp, double **xbounoldp, ITG *iexpl, double *plicon, ITG *nplicon, double *plkcon, ITG *nplkcon, double **xstatep, ITG *npmat_, char *matname, ITG *mi, ITG *ncmat_, ITG *nstate_, double **enerp, char *jobnamec, double *ttime, char *set, ITG *nset, ITG *istartset, ITG *iendset, ITG **ialsetp, ITG *nprint, char *prlab, char *prset, ITG *nener, double *trab, ITG **inotrp, ITG *ntrans, double **fmpcp, char *cbody, ITG *ibody, double *xbody, ITG *nbody, double *xbodyold, ITG *istep, ITG *isolver, ITG *jq, char *output, ITG *mcs, ITG *nkon, ITG *mpcend, ITG *ics, double *cs, ITG *ntie, char *tieset, ITG *idrct, ITG *jmax, double *ctrl, ITG *itpamp, double *tietol, ITG *nalset, ITG *ikforc, ITG *ilforc, double *thicke, ITG *nslavs, ITG *nmat, char *typeboun, ITG *ielprop, double *prop)
Definition: dyna.c:36
void complexfreq(double **cop, ITG *nk, ITG **konp, ITG **ipkonp, char **lakonp, ITG *ne, ITG **nodebounp, ITG **ndirbounp, double **xbounp, ITG *nboun, ITG **ipompcp, ITG **nodempcp, double **coefmpcp, char **labmpcp, ITG *nmpc, ITG *nodeforc, ITG *ndirforc, double *xforc, ITG *nforc, ITG *nelemload, char *sideload, double *xload, ITG *nload, ITG **nactdofp, ITG *neq, ITG *nzl, ITG *icol, ITG *irow, ITG *nmethod, ITG **ikmpcp, ITG **ilmpcp, ITG **ikbounp, ITG **ilbounp, double *elcon, ITG *nelcon, double *rhcon, ITG *nrhcon, double *cocon, ITG *ncocon, double *alcon, ITG *nalcon, double *alzero, ITG **ielmatp, ITG **ielorienp, ITG *norien, double *orab, ITG *ntmat_, double **t0p, double **t1p, ITG *ithermal, double *prestr, ITG *iprestr, double **voldp, ITG *iperturb, double **stip, ITG *nzs, double *timepar, double *xmodal, double **veoldp, char *amname, double *amta, ITG *namta, ITG *nam, ITG *iamforc, ITG *iamload, ITG **iamt1p, ITG *jout, ITG *kode, char *filab, double **emep, double *xforcold, double *xloadold, double **t1oldp, ITG **iambounp, double **xbounoldp, ITG *iexpl, double *plicon, ITG *nplicon, double *plkcon, ITG *nplkcon, double *xstate, ITG *npmat_, char *matname, ITG *mi, ITG *ncmat_, ITG *nstate_, double **enerp, char *jobnamec, double *ttime, char *set, ITG *nset, ITG *istartset, ITG *iendset, ITG **ialsetp, ITG *nprint, char *prlab, char *prset, ITG *nener, double *trab, ITG **inotrp, ITG *ntrans, double **fmpcp, char *cbody, ITG *ibody, double *xbody, ITG *nbody, double *xbodyold, ITG *istep, ITG *isolver, ITG *jq, char *output, ITG *mcs, ITG *nkon, ITG *mpcend, ITG *ics, double *cs, ITG *ntie, char *tieset, ITG *idrct, ITG *jmax, double *ctrl, ITG *itpamp, double *tietol, ITG *nalset, ITG *ikforc, ITG *ilforc, double *thicke, char *jobnamef, ITG *mei, ITG *nmat, ITG *ielprop, double *prop)
Definition: complexfreq.c:36
void nonlingeo(double **co, ITG *nk, ITG **konp, ITG **ipkonp, char **lakonp, ITG *ne, ITG *nodeboun, ITG *ndirboun, double *xboun, ITG *nboun, ITG **ipompcp, ITG **nodempcp, double **coefmpcp, char **labmpcp, 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 **ikmpcp, ITG **ilmpcp, ITG *ikboun, ITG *ilboun, double *elcon, ITG *nelcon, double *rhcon, ITG *nrhcon, double *alcon, ITG *nalcon, double *alzero, ITG **ielmatp, ITG **ielorienp, 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, ITG *idrct, ITG *jmax, ITG *jout, double *timepar, double *eme, double *xbounold, double *xforcold, double *xloadold, double *veold, double *accold, char *amname, double *amta, ITG *namta, ITG *nam, ITG *iamforc, ITG *iamload, ITG *iamt1, double *alpha, ITG *iexpl, ITG *iamboun, double *plicon, ITG *nplicon, double *plkcon, ITG *nplkcon, double **xstatep, ITG *npmat_, ITG *istep, double *ttime, char *matname, double *qaold, ITG *mi, ITG *isolver, ITG *ncmat_, ITG *nstate_, ITG *iumat, double *cs, ITG *mcs, ITG *nkon, double **ener, ITG *mpcinfo, char *output, double *shcon, ITG *nshcon, double *cocon, ITG *ncocon, double *physcon, ITG *nflow, double *ctrl, char *set, ITG *nset, ITG *istartset, ITG *iendset, ITG *ialset, ITG *nprint, char *prlab, char *prset, ITG *nener, ITG *ikforc, ITG *ilforc, double *trab, ITG *inotr, ITG *ntrans, double **fmpcp, char *cbody, ITG *ibody, double *xbody, ITG *nbody, double *xbodyold, ITG *ielprop, double *prop, ITG *ntie, char *tieset, ITG *itpamp, ITG *iviewfile, char *jobnamec, double *tietol, ITG *nslavs, double *thicke, ITG *ics, ITG *nintpoint, ITG *mortar, ITG *ifacecount, char *typeboun, ITG **islavsurfp, double **pslavsurfp, double **clearinip, ITG *nmat, double *xmodal, ITG *iaxial, ITG *inext)
Definition: nonlingeo.c:36
#define RENEW(a, b, c)
Definition: CalculiX.h:40
#define SFREE(a)
Definition: CalculiX.h:41
void writeheading(char *jobnamec, char *heading, ITG *nheading)
Definition: writeheading.c:24
subroutine allocation(nload, nforc, nboun, nk, ne, nmpc, nset, nalset, nmat, ntmat, npmat, norien, nam, nprint, mi, ntrans, set, meminset, rmeminset, ncs, namtot, ncmat, memmpc, ne1d, ne2d, nflow, jobnamec, irstrt, ithermal, nener, nstate, irestartstep, inpc, ipoinp, inp, ntie, nbody, nprop, ipoinpc, nevdamp, npt, nslavs, nkon, mcs, mortar, ifacecount, nintpoint, infree, nheading)
Definition: allocation.f:19
void mastructcs(ITG *nk, ITG *kon, ITG *ipkon, char *lakon, ITG *ne, ITG *nodeboun, ITG *ndirboun, ITG *nboun, ITG *ipompc, ITG *nodempc, ITG *nmpc, ITG *nactdof, ITG *icol, ITG *jq, ITG **mast1p, ITG **irowp, ITG *isolver, ITG *neq, ITG *ikmpc, ITG *ilmpc, ITG *ipointer, ITG *nzs, ITG *nmethod, ITG *ics, double *cs, char *labmpc, ITG *mcs, ITG *mi, ITG *mortar)
Definition: mastructcs.c:27
subroutine calinput(co, nk, kon, ipkon, lakon, nkon, ne, nodeboun, ndirboun, xboun, nboun, ipompc, nodempc, coefmpc, nmpc, nmpc_, nodeforc, ndirforc, xforc, nforc, nforc_, nelemload, sideload, xload, nload, nload_, nprint, prlab, prset, mpcfree, nboun_, mei, set, istartset, iendset, ialset, nset, nalset, elcon, nelcon, rhcon, nrhcon, alcon, nalcon, alzero, t0, t1, matname, ielmat, orname, orab, ielorien, amname, amta, namta, nam, nmethod, iamforc, iamload, iamt1, ithermal, iperturb, istat, istep, nmat, ntmat_, norien, prestr, iprestr, isolver, fei, veold, timepar, xmodal, filab, jout, nlabel, idrct, jmax, iexpl, alpha, iamboun, plicon, nplicon, plkcon, nplkcon, iplas, npmat_, mi, nk_, trab, inotr, ntrans, ikboun, ilboun, ikmpc, ilmpc, ics, dcs, ncs_, namtot_, cs, nstate_, ncmat_, iumat, mcs, labmpc, iponor, xnor, knor, thickn, thicke, ikforc, ilforc, offset, iponoel, inoel, rig, infree, nshcon, shcon, cocon, ncocon, physcon, nflow, ctrl, maxlenmpc, ne1d, ne2d, nener, vold, nodebounold, ndirbounold, xbounold, xforcold, xloadold, t1old, eme, sti, ener, xstate, jobnamec, irstrt, ttime, qaold, output, typeboun, inpc, ipoinp, inp, tieset, tietol, ntie, fmpc, cbody, ibody, xbody, nbody, nbody_, xbodyold, nam_, ielprop, nprop, nprop_, prop, itpamp, iviewfile, ipoinpc, cfd, nslavs, t0g, t1g, network, cyclicsymmetry, idefforc, idefload, idefbody, mortar, ifacecount, islavsurf, pslavsurf, clearini, heading, iaxial)
Definition: calinput.f:19
void arpackcs(double *co, ITG *nk, ITG **konp, ITG **ipkonp, char **lakonp, 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 *icol, 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 **ielmatp, ITG **ielorienp, 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, ITG *mei, double *fei, char *filab, double *eme, ITG *iexpl, double *plicon, ITG *nplicon, double *plkcon, ITG *nplkcon, double **xstatep, ITG *npmat_, char *matname, ITG *mi, ITG *ics, double *cs, ITG *mpcend, ITG *ncmat_, ITG *nstate_, ITG *mcs, ITG *nkon, double **enerp, char *jobnamec, char *output, char *set, ITG *nset, ITG *istartset, ITG *iendset, ITG *ialset, ITG *nprint, char *prlab, char *prset, ITG *nener, ITG *isolver, double *trab, ITG *inotr, ITG *ntrans, double *ttime, double *fmpc, char *cbody, ITG *ibody, double *xbody, ITG *nbody, ITG *nevtot, double *thicke, ITG *nslavs, double *tietol, ITG *mpcinfo, ITG *ntie, ITG *istep, char *tieset, ITG *nintpoint, ITG *mortar, ITG *ifacecount, ITG **islavsurfp, double **pslavsurfp, double **clearinip, ITG *nmat, char *typeboun, ITG *ielprop, double *prop)
void steadystate(double **co, ITG *nk, ITG **kon, ITG **ipkon, char **lakon, ITG *ne, ITG **nodeboun, ITG **ndirboun, double **xboun, ITG *nboun, ITG **ipompcp, ITG **nodempcp, double **coefmpcp, char **labmpcp, ITG *nmpc, ITG *nodeforc, ITG *ndirforc, double *xforc, ITG *nforc, ITG *nelemload, char *sideload, double *xload, ITG *nload, ITG **nactdof, ITG *neq, ITG *nzl, ITG *icol, ITG *irow, ITG *nmethod, ITG **ikmpcp, ITG **ilmpcp, ITG **ikboun, ITG **ilboun, double *elcon, ITG *nelcon, double *rhcon, ITG *nrhcon, double *cocon, ITG *ncocon, 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, double **voldp, ITG *iperturb, double *sti, ITG *nzs, double *timepar, double *xmodal, double **veoldp, char *amname, double *amta, ITG *namta, ITG *nam, ITG *iamforc, ITG *iamload, ITG **iamt1, ITG *jout, ITG *kode, char *filab, double **emep, double *xforcold, double *xloadold, double **t1old, ITG **iamboun, double **xbounold, ITG *iexpl, double *plicon, ITG *nplicon, double *plkcon, ITG *nplkcon, double *xstate, ITG *npmat_, char *matname, ITG *mi, ITG *ncmat_, ITG *nstate_, double **enerp, char *jobnamec, double *ttime, 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 **fmpcp, char *cbody, ITG *ibody, double *xbody, ITG *nbody, double *xbodyold, ITG *istep, ITG *isolver, ITG *jq, char *output, ITG *mcs, ITG *nkon, ITG *ics, double *cs, ITG *mpcend, double *ctrl, ITG *ikforc, ITG *ilforc, double *thicke, ITG *nmat, char *typeboun, ITG *ielprop, double *prop)
Definition: steadystate.c:37
void arpack(double *co, ITG *nk, ITG **konp, ITG **ipkonp, char **lakonp, 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 *icol, 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 *shcon, ITG *nshcon, double *cocon, ITG *ncocon, double *alcon, ITG *nalcon, double *alzero, ITG **ielmatp, ITG **ielorienp, 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, ITG *mei, double *fei, char *filab, double *eme, ITG *iexpl, double *plicon, ITG *nplicon, double *plkcon, ITG *nplkcon, double **xstatep, ITG *npmat_, char *matname, ITG *mi, ITG *ncmat_, ITG *nstate_, double **enerp, char *jobnamec, char *output, char *set, ITG *nset, ITG *istartset, ITG *iendset, ITG *ialset, ITG *nprint, char *prlab, char *prset, ITG *nener, ITG *isolver, double *trab, ITG *inotr, ITG *ntrans, double *ttime, double *fmpc, char *cbody, ITG *ibody, double *xbody, ITG *nbody, double *thicke, ITG *nslavs, double *tietol, ITG *nkon, ITG *mpcinfo, ITG *ntie, ITG *istep, ITG *mcs, ITG *ics, char *tieset, double *cs, ITG *nintpoint, ITG *mortar, ITG *ifacecount, ITG **islavsurfp, double **pslavsurfp, double **clearinip, ITG *nmat, char *typeboun, ITG *ielprop, double *prop)
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: linstatic.c:35
#define ITG
Definition: CalculiX.h:51
#define NNEW(a, b, c)
Definition: CalculiX.h:39
subroutine spcmatch(xboun, nodeboun, ndirboun, nboun, xbounold, nodebounold, ndirbounold, nbounold, ikboun, ilboun, vold, reorder, nreorder, mi)
Definition: spcmatch.f:19
subroutine closefile()
Definition: closefile.f:19
static ITG * nload1
Definition: results.c:31
subroutine genadvecelem(inodesd, ipkon, ne, lakon, kon, nload, sideload, nelemload, nkon)
Definition: genadvecelem.f:19
Hosted by OpenAircraft.com, (Michigan UAV, LLC)