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

Go to the source code of this file.

Functions

void compfluid (double **cop, ITG *nk, ITG **ipkonfp, ITG **konp, char **lakonfp, char **sidefacep, ITG *ifreestream, ITG *nfreestream, ITG *isolidsurf, ITG *neighsolidsurf, ITG *nsolidsurf, ITG **iponoelp, ITG **inoelp, ITG *nshcon, double *shcon, ITG *nrhcon, double *rhcon, double **voldp, ITG *ntmat_, ITG *nodeboun, ITG *ndirboun, ITG *nboun, ITG **ipompcp, ITG **nodempcp, ITG *nmpc, ITG **ikmpcp, ITG **ilmpcp, ITG *ithermal, ITG *ikboun, ITG *ilboun, ITG *iturbulent, ITG *isolver, ITG *iexpl, double *vcontu, double *ttime, double *time, double *dtime, ITG *nodeforc, ITG *ndirforc, double *xforc, ITG *nforc, ITG *nelemload, char *sideload, double *xload, ITG *nload, double *xbody, ITG *ipobody, ITG *nbody, ITG *ielmatf, char *matname, ITG *mi, ITG *ncmat_, double *physcon, ITG *istep, ITG *iinc, ITG *ibody, double *xloadold, double *xboun, double **coefmpcp, ITG *nmethod, double *xforcold, double *xforcact, ITG *iamforc, ITG *iamload, double *xbodyold, double *xbodyact, double *t1old, double *t1, double *t1act, ITG *iamt1, double *amta, ITG *namta, ITG *nam, double *ampli, double *xbounold, double *xbounact, ITG *iamboun, ITG *itg, ITG *ntg, char *amname, double *t0, ITG **nelemfacep, ITG *nface, double *cocon, ITG *ncocon, double *xloadact, double *tper, ITG *jmax, ITG *jout, char *set, ITG *nset, ITG *istartset, ITG *iendset, ITG *ialset, char *prset, char *prlab, ITG *nprint, double *trab, ITG *inotr, ITG *ntrans, char *filab, char **labmpcp, double *sti, ITG *norien, double *orab, char *jobnamef, char *tieset, ITG *ntie, ITG *mcs, ITG *ics, double *cs, ITG *nkon, ITG *mpcfree, ITG *memmpc_, double **fmpcp, ITG *nef, ITG **inomatp, double *qfx, ITG *neifa, ITG *neiel, ITG *ielfa, ITG *ifaext, double *vfa, double *vel, ITG *ipnei, ITG *nflnei, ITG *nfaext, char *typeboun, ITG *neij, double *tincf, ITG *nactdoh, ITG *nactdohinv, ITG *ielorienf)
 

Variables

static ITG num_cpus
 

Function Documentation

void compfluid ( double **  cop,
ITG nk,
ITG **  ipkonfp,
ITG **  konp,
char **  lakonfp,
char **  sidefacep,
ITG ifreestream,
ITG nfreestream,
ITG isolidsurf,
ITG neighsolidsurf,
ITG nsolidsurf,
ITG **  iponoelp,
ITG **  inoelp,
ITG nshcon,
double *  shcon,
ITG nrhcon,
double *  rhcon,
double **  voldp,
ITG ntmat_,
ITG nodeboun,
ITG ndirboun,
ITG nboun,
ITG **  ipompcp,
ITG **  nodempcp,
ITG nmpc,
ITG **  ikmpcp,
ITG **  ilmpcp,
ITG ithermal,
ITG ikboun,
ITG ilboun,
ITG iturbulent,
ITG isolver,
ITG iexpl,
double *  vcontu,
double *  ttime,
double *  time,
double *  dtime,
ITG nodeforc,
ITG ndirforc,
double *  xforc,
ITG nforc,
ITG nelemload,
char *  sideload,
double *  xload,
ITG nload,
double *  xbody,
ITG ipobody,
ITG nbody,
ITG ielmatf,
char *  matname,
ITG mi,
ITG ncmat_,
double *  physcon,
ITG istep,
ITG iinc,
ITG ibody,
double *  xloadold,
double *  xboun,
double **  coefmpcp,
ITG nmethod,
double *  xforcold,
double *  xforcact,
ITG iamforc,
ITG iamload,
double *  xbodyold,
double *  xbodyact,
double *  t1old,
double *  t1,
double *  t1act,
ITG iamt1,
double *  amta,
ITG namta,
ITG nam,
double *  ampli,
double *  xbounold,
double *  xbounact,
ITG iamboun,
ITG itg,
ITG ntg,
char *  amname,
double *  t0,
ITG **  nelemfacep,
ITG nface,
double *  cocon,
ITG ncocon,
double *  xloadact,
double *  tper,
ITG jmax,
ITG jout,
char *  set,
ITG nset,
ITG istartset,
ITG iendset,
ITG ialset,
char *  prset,
char *  prlab,
ITG nprint,
double *  trab,
ITG inotr,
ITG ntrans,
char *  filab,
char **  labmpcp,
double *  sti,
ITG norien,
double *  orab,
char *  jobnamef,
char *  tieset,
ITG ntie,
ITG mcs,
ITG ics,
double *  cs,
ITG nkon,
ITG mpcfree,
ITG memmpc_,
double **  fmpcp,
ITG nef,
ITG **  inomatp,
double *  qfx,
ITG neifa,
ITG neiel,
ITG ielfa,
ITG ifaext,
double *  vfa,
double *  vel,
ITG ipnei,
ITG nflnei,
ITG nfaext,
char *  typeboun,
ITG neij,
double *  tincf,
ITG nactdoh,
ITG nactdohinv,
ITG ielorienf 
)

Definition at line 39 of file compfluid.c.

67  {
68 
69  /* main computational fluid dynamics routine */
70 
71  char cflag[1],*labmpc=NULL,*lakonf=NULL,*sideface=NULL;
72 
73  char matvec[7]="MATVEC",msolve[7]="MSOLVE";
74 
75  ITG *ipointer=NULL,*mast1=NULL,*irow=NULL,*icol=NULL,*jq=NULL,
76  nzs=20000000,neq,kode,compressible,*ifabou=NULL,*ja=NULL,
77  *nodempc=NULL,*ipompc=NULL,*ikmpc=NULL,*ilmpc=NULL,nfabou,im,
78  *ipkonf=NULL,*kon=NULL,*nelemface=NULL,*inoel=NULL,last=0,
79  *iponoel=NULL,*inomat=NULL,ithermalref,*integerglob=NULL,iit,
80  iconvergence=0,symmetryflag,inputformat,i,*inum=NULL,iitf,ifreefa,
81  *iponofa=NULL,*inofa=NULL,is,ie,*ia=NULL,nstate_,*ielpropf=NULL,
82  icent=0,isti=0,iqfx=0,nfield,ndim,iorienglob,force=0,icfdout=1;
83 
84  ITG nelt,isym,itol,itmax,iunit,lrgw,*igwk=NULL,ligw,ierr,*iwork=NULL,iter,
85  nsave,lenw,leniw;
86 
87  double *coefmpc=NULL,*fmpc=NULL,*umfa=NULL,reltime,*doubleglob=NULL,
88  *co=NULL,*vold=NULL,*coel=NULL,*cosa=NULL,*gradvel=NULL,*gradvfa=NULL,
89  *xxn=NULL,*xxi=NULL,*xle=NULL,*xlen=NULL,*xlet=NULL,timef,dtimef,
90  *cofa=NULL,*area=NULL,*xrlfa=NULL,reltimef,ttimef,*hcfa=NULL,*cpel=NULL,
91  *au=NULL,*ad=NULL,*b=NULL,*volume=NULL,*body=NULL,sigma=0.,betam,
92  *adb=NULL,*aub=NULL,*advfa=NULL,*ap=NULL,*bp=NULL,*xxj=NULL,
93  *v=NULL,*velo=NULL,*veloo=NULL,*gammat=NULL,*cosb=NULL,dmin,tincfguess,
94  *hel=NULL,*hfa=NULL,*auv=NULL,*adv=NULL,*bv=NULL,*sel=NULL,*gamma=NULL,
95  *gradtfa=NULL,*gradtel=NULL,*umel=NULL,*cpfa=NULL,*gradpel=NULL,
96  *fn=NULL,*eei=NULL,*xstate=NULL,*ener=NULL,*thicke=NULL,*eme=NULL,
97  ptimef,*stn=NULL,*qfn=NULL,*hcel=NULL,*aua=NULL,a1,a2,a3,beta,
98  *prop=NULL;
99 
100  double tol,*rgwk=NULL,err,*sb=NULL,*sx=NULL,*rwork=NULL;
101 
102  nodempc=*nodempcp;ipompc=*ipompcp;ikmpc=*ikmpcp;ilmpc=*ilmpcp;
103  coefmpc=*coefmpcp;labmpc=*labmpcp;fmpc=*fmpcp;co=*cop;
104  ipkonf=*ipkonfp;lakonf=*lakonfp;kon=*konp;
105  nelemface=*nelemfacep;sideface=*sidefacep;inoel=*inoelp;
106  iponoel=*iponoelp;vold=*voldp;inomat=*inomatp;
107 
108 #ifdef SGI
109  ITG token;
110 #endif
111 
112  /* relative time at the end of the mechanical increment */
113 
114  reltime=(*time)/(*tper);
115 
116  /* open frd-file for fluids */
117 
118  FORTRAN(openfilefluid,(jobnamef));
119 
120  /* variables for multithreading procedure */
121 
122  ITG sys_cpus;
123  char *env,*envloc,*envsys;
124 
125  num_cpus = 0;
126  sys_cpus=0;
127 
128  /* explicit user declaration prevails */
129 
130  envsys=getenv("NUMBER_OF_CPUS");
131  if(envsys){
132  sys_cpus=atoi(envsys);
133  if(sys_cpus<0) sys_cpus=0;
134  }
135 
136  /* automatic detection of available number of processors */
137 
138  if(sys_cpus==0){
139  sys_cpus = getSystemCPUs();
140  if(sys_cpus<1) sys_cpus=1;
141  }
142 
143  /* local declaration prevails, if strictly positive */
144 
145  envloc = getenv("CCX_NPROC_CFD");
146  if(envloc){
147  num_cpus=atoi(envloc);
148  if(num_cpus<0){
149  num_cpus=0;
150  }else if(num_cpus>sys_cpus){
151  num_cpus=sys_cpus;
152  }
153  }
154 
155  /* else global declaration, if any, applies */
156 
157  env = getenv("OMP_NUM_THREADS");
158  if(num_cpus==0){
159  if (env)
160  num_cpus = atoi(env);
161  if (num_cpus < 1) {
162  num_cpus=1;
163  }else if(num_cpus>sys_cpus){
164  num_cpus=sys_cpus;
165  }
166  }
167 
168 // next line is to be inserted in a similar way for all other paralell parts
169 
170  if(*nef<num_cpus) num_cpus=*nef;
171 
172  printf(" Using up to %" ITGFORMAT " cpu(s) for CFD.\n", num_cpus);
173 
174  pthread_t tid[num_cpus];
175 
176 
177  kode=0;
178 
179  /* *iexpl==0: structure:implicit, fluid:incompressible
180  *iexpl==1: structure:implicit, fluid:compressible
181  *iexpl==2: structure:explicit, fluid:incompressible
182  *iexpl==3: structure:explicit, fluid:compressible */
183 
184  if((*iexpl==1)||(*iexpl==3)){
185  compressible=1;
186  }else{
187  compressible=0;
188  }
189 
190  /* if initial conditions are specified for the temperature,
191  it is assumed that the temperature is an unknown */
192 
193  ithermalref=*ithermal;
194  if(*ithermal==1){
195  *ithermal=2;
196  }
197 
198  /* determining the matrix structure */
199 
200  NNEW(ipointer,ITG,3**nk);
201  NNEW(mast1,ITG,nzs);
202  NNEW(irow,ITG,nzs);
203  NNEW(ia,ITG,nzs);
204  NNEW(icol,ITG,*nef);
205  NNEW(jq,ITG,*nef+1);
206  NNEW(ja,ITG,*nef+1);
207 // NNEW(nactdoh,ITG,*nef);
208 
209  mastructf(nk,kon,ipkonf,lakonf,nef,icol,jq,&mast1,&irow,
210  isolver,&neq,ipointer,&nzs,ipnei,neiel,mi);
211 
212 // printf("Unterschied start\n");
213 // for(i=0;i<*ne;i++){if(i+1!=nactdoh[i]){printf("Unterschied i=%d,nactdoh[i]=%d\n",i+1,nactdoh[i]);}}
214 // printf("Unterschied end\n");
215 
216  SFREE(ipointer);SFREE(mast1);
217 
218  /* calculation geometric data */
219 
220  NNEW(coel,double,3**nef);
221  NNEW(volume,double,*nef);
222  NNEW(cosa,double,*nflnei);
223  NNEW(cosb,double,*nflnei);
224  NNEW(xxn,double,3**nflnei);
225  NNEW(xxi,double,3**nflnei);
226  NNEW(xxj,double,3**nflnei);
227  NNEW(xle,double,*nflnei);
228  NNEW(xlen,double,*nflnei);
229  NNEW(xlet,double,*nflnei);
230  NNEW(cofa,double,3**nface);
231  NNEW(area,double,*nface);
232  NNEW(xrlfa,double,3**nface);
233 
234  FORTRAN(initialcfd,(nef,ipkonf,kon,lakonf,co,coel,cofa,nface,
235  ielfa,area,ipnei,neiel,xxn,xxi,xle,xlen,xlet,xrlfa,cosa,
236  volume,neifa,xxj,cosb,vel,&dmin));
237 
238  /* storing pointers to the boundary conditions in ielfa */
239 
240  NNEW(ifabou,ITG,7**nfaext);
241  FORTRAN(applyboun,(ifaext,nfaext,ielfa,ikboun,ilboun,
242  nboun,typeboun,nelemload,nload,sideload,isolidsurf,nsolidsurf,
243  ifabou,&nfabou,nface,nodeboun,ndirboun,ikmpc,ilmpc,labmpc,nmpc,
244  nactdohinv));
245  RENEW(ifabou,ITG,nfabou);
246 
247  /* catalogueing the nodes for output purposes (interpolation at
248  the nodes */
249 
250  NNEW(iponofa,ITG,*nk);
251  NNEW(inofa,ITG,2**nface*4);
252 
253  FORTRAN(cataloguenodes,(iponofa,inofa,&ifreefa,ielfa,ifabou,ipkonf,
254  kon,lakonf,nface,nk));
255 
256  RENEW(inofa,ITG,2*ifreefa);
257 
258  /* material properties for athermal calculations
259  = calculation for which no initial thermal conditions
260  were defined */
261 
262  NNEW(umfa,double,*nface);
263  NNEW(umel,double,*nef);
264 
265  /* calculating the density at the element centers */
266 
267  FORTRAN(calcrhoel,(nef,vel,rhcon,nrhcon,ielmatf,ntmat_,
268  ithermal,mi));
269 
270  /* calculating the density at the face centers */
271 
272  FORTRAN(calcrhofa,(nface,vfa,rhcon,nrhcon,ielmatf,ntmat_,
273  ithermal,mi,ielfa));
274 
275  /* calculating the dynamic viscosity at the face centers */
276 
277  FORTRAN(calcumfa,(nface,vfa,shcon,nshcon,ielmatf,ntmat_,
278  ithermal,mi,ielfa,umfa));
279 
280  /* calculating the dynamic viscosity at the element centers */
281 
282  FORTRAN(calcumel,(nef,vel,shcon,nshcon,ielmatf,ntmat_,
283  ithermal,mi,umel));
284 
285  if(*ithermal!=0){
286  NNEW(hcfa,double,*nface);
287  NNEW(cpel,double,*nef);
288  NNEW(cpfa,double,*nface);
289  }
290 
291 
292  if(*nbody>0) NNEW(body,double,4**nef);
293 
294  /* v is a auxiliary field: set to zero for the calls to
295  tempload */
296 
297  NNEW(v,double,5**nk);
298 
299  /* next section is for stationary calculations */
300 
301  if(*nmethod==1){
302 
303  /* boundary conditions at the end of the mechanical
304  increment */
305 
306  FORTRAN(tempload,(xforcold,xforc,xforcact,iamforc,nforc,
307  xloadold,xload,xloadact,iamload,nload,ibody,xbody,nbody,
308  xbodyold,xbodyact,t1old,t1,t1act,iamt1,nk,amta,
309  namta,nam,ampli,time,&reltime,ttime,dtime,ithermal,nmethod,
310  xbounold,xboun,xbounact,iamboun,nboun,
311  nodeboun,ndirboun,nodeforc,ndirforc,istep,iinc,
312  co,v,itg,ntg,amname,ikboun,ilboun,nelemload,sideload,mi,
313  ntrans,trab,inotr,vold,integerglob,doubleglob,tieset,istartset,
314  iendset,ialset,ntie,nmpc,ipompc,ikmpc,ilmpc,nodempc,coefmpc));
315 
316  /* body forces (gravity, centrifugal and Coriolis forces */
317 
318  if(*nbody>0){
319  FORTRAN(inicalcbody,(nef,body,ipobody,ibody,xbody,coel,vel,lakonf,
320  nactdohinv,&icent));
321  }
322  }
323 
324  /* extrapolating the velocity from the elements centers to the face
325  centers, thereby taking the boundary conditions into account */
326 
327  FORTRAN(extrapolate_vel,(nface,ielfa,xrlfa,vel,vfa,
328  ifabou,xbounact,ipnei,nef));
329 
330  /* applying MPC's to the faces */
331 
332  if(*nmpc>0){
333  is=1;ie=3;
334  FORTRAN(applympc,(nface,ielfa,&is,&ie,ifabou,ipompc,vfa,coefmpc,
335  nodempc,ipnei,neifa,labmpc,xbounact,nactdoh));
336  }
337 
338  /* extrapolation of the pressure at the element centers
339  to the face centers */
340 
341  FORTRAN(extrapolate_pel,(nface,ielfa,xrlfa,vel,vfa,ifabou,
342  xbounact,nef));
343 
344  /* applying MPC's to the faces */
345 
346  if(*nmpc>0){
347  is=4;ie=4;
348  FORTRAN(applympc,(nface,ielfa,&is,&ie,ifabou,ipompc,vfa,coefmpc,
349  nodempc,ipnei,neifa,labmpc,xbounact,nactdoh));
350  }
351 
352  /* extrapolation of the temperature at the element centers
353  to the face centers */
354 
355  if(*ithermal>0){
356  FORTRAN(extrapolate_tel,(nface,ielfa,xrlfa,vel,vfa,
357  ifabou,xbounact,ipnei,nef));
358 
359  /* applying MPC's to the faces */
360 
361  if(*nmpc>0){
362  is=0;ie=0;
363  FORTRAN(applympc,(nface,ielfa,&is,&ie,ifabou,ipompc,vfa,coefmpc,
364  nodempc,ipnei,neifa,labmpc,xbounact,nactdoh));
365  }
366  }
367 
368  /* calculating the maximum velocity (for the determination of the
369  time increment) */
370 
371  FORTRAN(calcguesstincf,(nface,&dmin,vfa,umfa,&tincfguess));
372 
373  /* start of the major loop */
374 
375  NNEW(gradvel,double,9**nef);
376  NNEW(gradvfa,double,9**nface);
377 
378  NNEW(advfa,double,*nface);
379  NNEW(hfa,double,3**nface);
380 
381  NNEW(ap,double,*nface);
382  NNEW(bp,double,*nface);
383 
384  NNEW(au,double,2*nzs+neq);
385  NNEW(aua,double,nzs+neq);
386 // NNEW(au,double,2*nzs);
387  NNEW(ad,double,neq);
388  NNEW(b,double,neq);
389 
390  NNEW(auv,double,2*nzs+neq);
391 // NNEW(auv,double,2*nzs);
392  NNEW(adv,double,neq);
393  NNEW(bv,double,3*neq);
394  NNEW(hel,double,3*neq);
395  NNEW(sel,double,3*neq);
396 
397  NNEW(rwork,double,neq);
398 
399  NNEW(gradpel,double,3**nef);
400 
401  if(*ithermal>0){
402  NNEW(gradtel,double,3**nef);
403  NNEW(gradtfa,double,3**nface);
404  }
405 
406  NNEW(inum,ITG,*nk);
407 
408  NNEW(velo,double,6**nef);
409  NNEW(veloo,double,6**nef);
410 
411  /* initializing velo and veloo */
412 
413  memcpy(&veloo[0],&vel[0],sizeof(double)*6**nef);
414  memcpy(&velo[0],&vel[0],sizeof(double)*6**nef);
415 
416  iit=0;
417 
418  a1=1.5;
419  a2=-2.;
420  a3=0.5;
421 
422  if(*tincf<=0.) *tincf=tincfguess;
423  printf("time increment for the CFD-calculations = %e\n\n",*tincf);
424 
425  ttimef=*ttime;
426  timef=*time-*dtime;
427  dtimef=*tincf;
428 
429  do{
430 
431  iit++;
432 
433  printf("iteration = %d\n",iit);
434 
435  timef+=dtimef;
436  if((*time<timef)&&(*nmethod==4)){
437  dtimef-=(timef-*time);
438  timef=*time;
439  last=1;
440  beta=dtimef/(*tincf);
441  a1=(2.+beta)/(1.+beta);
442  a2=-(1.+beta)/beta;
443  a3=1./(beta*(1.+beta));
444  }
445 
446  /* conditions for transient calculations */
447 
448  if(*nmethod==4){
449 
450  /* boundary conditions at end of fluid increment */
451 
452  FORTRAN(tempload,(xforcold,xforc,xforcact,iamforc,nforc,
453  xloadold,xload,xloadact,iamload,nload,ibody,xbody,nbody,
454  xbodyold,xbodyact,t1old,t1,t1act,iamt1,nk,amta,
455  namta,nam,ampli,&timef,&reltimef,&ttimef,&dtimef,ithermal,nmethod,
456  xbounold,xboun,xbounact,iamboun,nboun,
457  nodeboun,ndirboun,nodeforc,ndirforc,istep,iinc,
458  co,v,itg,ntg,amname,ikboun,ilboun,nelemload,sideload,mi,
459  ntrans,trab,inotr,vold,integerglob,doubleglob,tieset,istartset,
460  iendset,ialset,ntie,nmpc,ipompc,ikmpc,ilmpc,nodempc,coefmpc));
461 
462  /* body forces (gravity, centrifugal and Coriolis forces) */
463 
464  if(*nbody>0){
465  FORTRAN(calcbody,(nef,body,ipobody,ibody,xbody,coel,vel,lakonf,
466  nactdohinv));
467  }
468 
469  }else if(icent==1){
470 
471  /* body forces (gravity, centrifugal and Coriolis forces;
472  only if centrifugal forces are active => the ensuing
473  Coriolis forces depend on the actual velocity) */
474 
475  FORTRAN(calcbody,(nef,body,ipobody,ibody,xbody,coel,vel,lakonf,
476  nactdohinv));
477  }
478 
479  /* updating of the material properties */
480 
481  if(*ithermal>0){
482 
483  /* calculating the density at the element centers */
484 
485  FORTRAN(calcrhoel,(nef,vel,rhcon,nrhcon,ielmatf,ntmat_,
486  ithermal,mi));
487 
488  /* calculating the density at the face centers */
489 
490  FORTRAN(calcrhofa,(nface,vfa,rhcon,nrhcon,ielmatf,ntmat_,
491  ithermal,mi,ielfa));
492 
493  /* calculating the dynamic viscosity at the face centers */
494 
495  FORTRAN(calcumfa,(nface,vfa,shcon,nshcon,ielmatf,ntmat_,
496  ithermal,mi,ielfa,umfa));
497 
498  /* calculating the dynamic viscosity at the element centers */
499 
500  FORTRAN(calcumel,(nef,vel,shcon,nshcon,ielmatf,ntmat_,
501  ithermal,mi,umel));
502 
503  /* calculating the heat conduction at the face centers */
504 
505  FORTRAN(calchcfa,(nface,vfa,cocon,ncocon,ielmatf,ntmat_,
506  mi,ielfa,hcfa));
507 
508  /* calculating the specific heat at constant pressure at the
509  element centers (secant value) */
510 
511  FORTRAN(calccpel,(nef,vel,shcon,nshcon,ielmatf,ntmat_,
512  mi,cpel,physcon));
513 
514  /* calculating the specific heat at constant pressure at the
515  face centers (secant value) */
516 
517  FORTRAN(calccpfa,(nface,vfa,shcon,nshcon,ielmatf,ntmat_,
518  mi,ielfa,cpfa,physcon));
519 
520  }
521 
522  /* calculating the gradient of the velocity at the element
523  centers */
524 
525  DMEMSET(gradvel,0,9**nef,0.);
526  FORTRAN(calcgradvel,(nef,lakonf,ipnei,vfa,area,xxn,gradvel,neifa,
527  volume));
528 
529  /* extrapolating the gradient of the velocity from the element
530  centers to the face centers */
531 
532  FORTRAN(extrapolate_gradvel,(nface,ielfa,xrlfa,gradvel,gradvfa));
533 
534  /* calculate gamma (Ph.D. Thesis Jasak) */
535 
536  /* betam=0.1;
537  NNEW(gamma,double,*nface);
538  FORTRAN(calcgamma,(nface,ielfa,vel,gradvfa,gamma,xlet,xxn,xxj,
539  ipnei,&betam,nef));*/
540 
541  /* filling the lhs and rhs's for the balance of momentum
542  equations */
543 
544  DMEMSET(adv,0,neq,0.);
545  DMEMSET(auv,0,2*nzs,0.);
546  DMEMSET(bv,0,3*neq,0.);
547 
548  FORTRAN(mafillv,(nef,ipnei,neifa,neiel,vfa,xxn,area,
549  auv,adv,jq,irow,&nzs,bv,vel,cosa,umfa,xlet,xle,gradvfa,xxi,
550  body,volume,&compressible,ielfa,lakonf,ifabou,nbody,&neq,
551  &dtimef,velo,veloo,sel,xrlfa,gamma,xxj,nactdohinv,&a1,
552  &a2,&a3));
553 // SFREE(gamma);
554 
555  /* LU decomposition (asymmetric system) */
556 
557 /* inputformat=1;
558  symmetryflag=2;
559 
560  if(*isolver==0){
561 #ifdef SPOOLES
562  spooles_factor(adv,auv,adb,aub,&sigma,icol,irow,&neq,&nzs,
563  &symmetryflag,&inputformat,&nzs);
564 #else
565  printf("*ERROR in compfluid: the SPOOLES library is not linked\n\n");
566  FORTRAN(stop,());
567 #endif
568  }
569  else if(*isolver==7){
570 #ifdef PARDISO
571  pardiso_factor(adv,auv,adb,aub,&sigma,icol,irow,&neq,&nzs,
572  &symmetryflag,&inputformat,jq,&nzs);
573 #else
574  printf("*ERROR in compfluid: the PARDISO library is not linked\n\n");
575  FORTRAN(stop,());
576 #endif
577 }*/
578 
579  /* solving the system of equations (for x, y and z
580  separately */
581 
582 /* for(i=0;i<3;i++){
583 
584  if(*isolver==0){
585 #ifdef SPOOLES
586  spooles_solve(&bv[i*neq],&neq);
587 #endif
588  }
589  else if(*isolver==7){
590 #ifdef PARDISO
591  pardiso_solve(&bv[i*neq],&neq,&symmetryflag);
592 #endif
593  }
594  }*/
595 
596  /* free memory */
597 
598 /* if(*isolver==0){
599 #ifdef SPOOLES
600  spooles_cleanup();
601 #endif
602  }
603  else if(*isolver==7){
604 #ifdef PARDISO
605  pardiso_cleanup(&neq,&symmetryflag);
606 #endif
607 }*/
608 
609  memcpy(&auv[2*nzs],adv,sizeof(double)*neq);
610  nelt=2*nzs+neq;
611  isym=0;
612  itol=0;
613  tol=0.;
614  itmax=0;
615  iunit=0;
616  lrgw=131+16*neq;
617  NNEW(rgwk,double,lrgw);
618  NNEW(igwk,ITG,20);
619  igwk[0]=10;
620  igwk[1]=10;
621  igwk[2]=0;
622  igwk[3]=1;
623  igwk[4]=10;
624  ligw=20;
625  for(i=0;i<neq;i++){rwork[i]=1./adv[i];}
626  FORTRAN(predgmres,(&neq,&bv[0],&vel[neq],&nelt,irow,jq,auv,
627  &isym,&itol,&tol,&itmax,&iter,
628  &err,&ierr,&iunit,sb,sx,rgwk,&lrgw,igwk,
629  &ligw,rwork,iwork));
630 // printf(" err1=%e\n",err);
631 // memcpy(&bv[0],&vel[neq],sizeof(double)*neq);
632  FORTRAN(predgmres,(&neq,&bv[neq],&vel[2*neq],&nelt,irow,jq,auv,
633  &isym,&itol,&tol,&itmax,&iter,
634  &err,&ierr,&iunit,sb,sx,rgwk,&lrgw,igwk,
635  &ligw,rwork,iwork));
636 // printf(" err2=%e\n",err);
637 // memcpy(&bv[neq],&vel[2*neq],sizeof(double)*neq);
638  FORTRAN(predgmres,(&neq,&bv[2*neq],&vel[3*neq],&nelt,irow,jq,auv,
639  &isym,&itol,&tol,&itmax,&iter,
640  &err,&ierr,&iunit,sb,sx,rgwk,&lrgw,igwk,
641  &ligw,rwork,iwork));
642 // printf(" err3=%e\n",err);
643 // memcpy(&bv[2*neq],&vel[3*neq],sizeof(double)*neq);
644  SFREE(rgwk);SFREE(igwk);
645 
646  /* storing the solution into vel */
647 
648 // FORTRAN(calcvel,(ne,nactdoh,vel,bv,&neq,ne));
649 
650  /* calculating the pressure gradient at the element
651  centers */
652 
653  DMEMSET(gradpel,0,3**nef,0.);
654  FORTRAN(calcgradpel,(nef,lakonf,ipnei,vfa,area,xxn,gradpel,neifa,
655  volume));
656 
657  for(iitf=0;iitf<2;iitf++){
658 
659  memcpy(&hel[0],&sel[0],sizeof(double)*(3*neq));
660 
661  /* completing hel with the neighboring velocity contributions */
662 
663  FORTRAN(complete_hel,(&neq,&vel[neq],hel,adv,auv,jq,irow,&nzs));
664 // FORTRAN(complete_hel,(&neq,bv,hel,adv,auv,jq,irow,&nzs));
665 
666  /* generating ad and h at the face centers (advfa and hfa) */
667 
668  FORTRAN(extrapolate_ad_h,(nface,ielfa,xrlfa,adv,advfa,hel,hfa));
669 
670  /* calculating the lhs and rhs of the equation system to determine
671  p (balance of mass) */
672 
673  DMEMSET(b,0,neq,0.);
674 
675  if(iitf==0){
676 
677  /* first iteration: calculating both lhs and rhs */
678 
679  DMEMSET(ad,0,neq,0.);
680  DMEMSET(au,0,nzs,0.);
681 
682  FORTRAN(mafillp,(nef,lakonf,ipnei,neifa,neiel,vfa,area,
683  advfa,xlet,cosa,volume,au,ad,jq,irow,ap,
684  ielfa,ifabou,xle,b,xxn,&compressible,&neq,
685  &nzs,hfa,gradpel,bp,xxi,neij,xlen,cosb));
686 
687  FORTRAN(convert2slapcol,(au,ad,irow,ia,jq,ja,&nzs,&neq,aua));
688 
689  /* LU decomposition of the p system (symmetric system) */
690 
691 /* inputformat=0;
692  symmetryflag=0;
693 
694  if(*isolver==0){
695 #ifdef SPOOLES
696  spooles_factor(ad,au,adb,aub,&sigma,icol,irow,&neq,&nzs,
697  &symmetryflag,&inputformat,&nzs);
698 #else
699  printf("*ERROR in compfluid: the SPOOLES library is not linked\n\n");
700  FORTRAN(stop,());
701 #endif
702  }
703  else if(*isolver==7){
704 #ifdef PARDISO
705  pardiso_factor(ad,au,adb,aub,&sigma,icol,irow,&neq,&nzs,
706  &symmetryflag,&inputformat,jq,&nzs);
707 #else
708  printf("*ERROR in compfluid: the PARDISO library is not linked\n\n");
709  FORTRAN(stop,());
710 #endif
711 }*/
712  }else{
713 
714  /* calculating the pressure gradient at the element
715  centers */
716 
717  DMEMSET(gradpel,0,3**nef,0.);
718  FORTRAN(calcgradpel,(nef,lakonf,ipnei,vfa,area,xxn,gradpel,neifa,
719  volume));
720 
721  /* second, third.. iteration: calculate the rhs only */
722 
723  FORTRAN(rhsp,(nef,lakonf,nactdoh,ipnei,neifa,neiel,vfa,area,
724  advfa,xlet,cosa,volume,au,ad,jq,irow,ap,ielfa,ifabou,xle,
725  b,xxn,&compressible,&neq,&nzs,hfa,bp,neij,xxi,
726  gradpel,xlen));
727 
728  }
729 
730  /* solving the system of equations for p */
731 
732 /* if(*isolver==0){
733 #ifdef SPOOLES
734  spooles_solve(b,&neq);
735 #endif
736  }
737  else if(*isolver==7){
738 #ifdef PARDISO
739  pardiso_solve(b,&neq,&symmetryflag);
740 #endif
741 }*/
742 
743  /* dslugm; attention: convert2slapcol changes au! */
744 
745 // FORTRAN(convert2slapcol,(au,ad,irow,ia,jq,ja,&nzs,&neq,aua));
746 
747  nelt=nzs+neq;
748  isym=1;
749  nsave=10;
750  itol=0;
751  tol=0.;
752  itmax=110;
753  iunit=0;
754  lenw=131+17*neq+2*nelt;
755  NNEW(rgwk,double,lenw);
756  leniw=32+4*neq+2*nelt;
757  NNEW(igwk,ITG,leniw);
758 
759  FORTRAN(dslugm,(&neq,&b[0],&vel[4*neq],&nelt,ia,ja,aua,
760  &isym,&nsave,&itol,&tol,&itmax,&iter,
761  &err,&ierr,&iunit,rgwk,&lenw,igwk,&leniw));
762  SFREE(rgwk);SFREE(igwk);
763 
764  if(ierr!=0){
765  printf(" error message from dslugm: ierr = %d\n",ierr);
766  printf(" err=%e\n",err);
767  }
768 
769  /* storing the solution p into vel(4,*) */
770 
771 // FORTRAN(calcpel,(nef,nactdoh,vel,b,nef));
772 
773  /* extrapolation of the pressure at the element centers
774  to the face centers */
775 
776  FORTRAN(extrapolate_pel,(nface,ielfa,xrlfa,vel,vfa,ifabou,
777  xbounact,nef));
778 
779  /* applying MPC's to the faces */
780 
781  if(*nmpc>0){
782  is=4;ie=4;
783  FORTRAN(applympc,(nface,ielfa,&is,&ie,ifabou,ipompc,vfa,coefmpc,
784  nodempc,ipnei,neifa,labmpc,xbounact,nactdoh));
785  }
786 
787  /* correction of the velocity at the element centers due
788  to the pressure change */
789 
790  FORTRAN(correctvel,(hel,adv,vfa,ipnei,area,&vel[neq],xxn,neifa,
791  lakonf,nef,&neq));
792 
793  }
794 
795  /* free memory */
796 
797  /* if(*isolver==0){
798 #ifdef SPOOLES
799  spooles_cleanup();
800 #endif
801  }
802  else if(*isolver==7){
803 #ifdef PARDISO
804  pardiso_cleanup(&neq,&symmetryflag);
805 #endif
806 }*/
807 
808  if(*ithermal>0){
809 
810  /* adding the velocity correction at the face centers
811  due to the balance of mass =>
812  the resulting mass flux is correct,
813  the face velocity vectors do not have to be correct
814  only needed if extra equations (temperature,
815  turbulence have to be solved) */
816 
817  FORTRAN(correctvfa,(nface,ielfa,area,vfa,ap,bp,xxn,
818  ifabou,ipnei,nef,neifa,hfa,vel,xbounact));
819 
820  /* calculating the temperature gradient at the element
821  centers */
822 
823  DMEMSET(gradtel,0,3**nef,0.);
824  FORTRAN(calcgradtel,(nef,lakonf,ipnei,vfa,area,xxn,gradtel,neifa,
825  volume));
826 
827  /* extrapolating the temperature gradient from the element
828  centers to the face centers */
829 
830  FORTRAN(extrapolate_gradtel,(nface,ielfa,xrlfa,gradtel,gradtfa));
831 
832  /* calculate gammat (Ph.D. Thesis Jasak) */
833 
834 // betam=0.1;
835 // NNEW(gammat,double,*nface);
836 // FORTRAN(calcgammat,(nface,ielfa,vel,gradtfa,gammat,xlet,xxn,xxj,
837 // ipnei,&betam,ne));
838 
839  /* calculating the lhs and rhs of the energy equation */
840 
841  DMEMSET(ad,0,neq,0.);
842  DMEMSET(au,0,2*nzs,0.);
843  DMEMSET(b,0,neq,0.);
844 
845  FORTRAN(mafillt,(nef,ipnei,neifa,neiel,vfa,xxn,area,
846  au,ad,jq,irow,&nzs,b,vel,umel,xlet,xle,gradtfa,xxi,
847  body,volume,&compressible,ielfa,lakonf,ifabou,nbody,&neq,
848  &dtimef,velo,veloo,cpfa,hcfa,cpel,gradvel,xload,gammat,
849  xrlfa,xxj,nactdohinv,&a1,&a2,&a3));
850 
851 // SFREE(gammat);
852 
853  /* solving the asymmetric system of equations */
854 
855  /* inputformat=1;
856  symmetryflag=2;
857 
858  if(*isolver==0){
859 #ifdef SPOOLES
860  spooles(ad,au,adb,aub,&sigma,b,icol,irow,&neq,&nzs,
861  &symmetryflag,&inputformat,&nzs);
862 #else
863  printf("*ERROR in compfluid: the SPOOLES library is not linked\n\n");
864  FORTRAN(stop,());
865 #endif
866  }
867  else if(*isolver==7){
868 #ifdef PARDISO
869  pardiso_main(ad,au,adb,aub,&sigma,b,icol,irow,&neq,&nzs,
870  &symmetryflag,&inputformat,jq,&nzs);
871 #else
872  printf("*ERROR in compfluid: the PARDISO library is not linked\n\n");
873  FORTRAN(stop,());
874 #endif
875 }*/
876 
877  /* free memory */
878 
879  /* if(*isolver==0){
880 #ifdef SPOOLES
881  spooles_cleanup();
882 #endif
883  }
884  else if(*isolver==7){
885 #ifdef PARDISO
886  pardiso_cleanup(&neq,&symmetryflag);
887 #endif
888 }*/
889 
890  memcpy(&au[2*nzs],ad,sizeof(double)*neq);
891  nelt=2*nzs+neq;
892  isym=0;
893  itol=0;
894  tol=0.;
895  itmax=0;
896  iunit=0;
897  lrgw=131+16*neq;
898  NNEW(rgwk,double,lrgw);
899  NNEW(igwk,ITG,20);
900  igwk[0]=10;
901  igwk[1]=10;
902  igwk[2]=0;
903  igwk[3]=1;
904  igwk[4]=10;
905  ligw=20;
906  for(i=0;i<neq;i++){rwork[i]=1./ad[i];}
907  FORTRAN(predgmres,(&neq,&b[0],&vel[0],&nelt,irow,jq,au,
908  &isym,&itol,&tol,&itmax,&iter,
909  &err,&ierr,&iunit,sb,sx,rgwk,&lrgw,igwk,
910  &ligw,rwork,iwork));
911  SFREE(rgwk);SFREE(igwk);
912  if(ierr>0){
913  printf("*WARNING in compfluid: error message from predgmres=%e\n",err);
914  }
915 
916  /* storing the solution t into vel(0,*) */
917 
918 // FORTRAN(calctel,(ne,nactdoh,vel,b,ne));
919 
920  /* extrapolation of the temperature at the element centers
921  to the face centers */
922 
923  FORTRAN(extrapolate_tel,(nface,ielfa,xrlfa,vel,vfa,
924  ifabou,xbounact,ipnei,nef));
925 
926  /* applying MPC's to the faces */
927 
928  if(*nmpc>0){
929  is=0;ie=0;
930  FORTRAN(applympc,(nface,ielfa,&is,&ie,ifabou,ipompc,vfa,coefmpc,
931  nodempc,ipnei,neifa,labmpc,xbounact,nactdoh));
932  }
933 
934  }
935 
936  /* storing the solution into vel */
937 
938 // FORTRAN(calcvel,(ne,nactdoh,vel,bv,&neq,ne));
939 
940  /* extrapolating the velocity from the elements centers to the face
941  centers, thereby taking the boundary conditions into account */
942 
943  FORTRAN(extrapolate_vel,(nface,ielfa,xrlfa,vel,vfa,
944  ifabou,xbounact,ipnei,nef));
945 
946  /* applying MPC's to the faces */
947 
948  if(*nmpc>0){
949  is=1;ie=3;
950  FORTRAN(applympc,(nface,ielfa,&is,&ie,ifabou,ipompc,vfa,coefmpc,
951  nodempc,ipnei,neifa,labmpc,xbounact,nactdoh));
952  }
953 
954 // FORTRAN(writevfa,(vfa,nface,nactdohinv,ielfa));
955 
956 
957  if(((iit/jout[1])*jout[1]==iit)||(iconvergence==1)||
958  (iit==jmax[1])){
959 
960  /* calculating the stress and the heat flow at the
961  integration points, if requested */
962 
963  if((strcmp1(&filab[3306],"SF ")==0)||
964  (strcmp1(&filab[3480],"SVF ")==0))isti=1;
965  if(strcmp1(&filab[3393],"HFLF")==0)iqfx=1;
966  for(i=0;i<*nprint;i++){
967  if(strcmp1(&prlab[6*i],"SVF")==0) isti=1;
968  if(strcmp1(&prlab[6*i],"HFLF")==0)iqfx=1;
969  }
970 
971  /* calculating the heat conduction at the element centers */
972 
973  if(iqfx==1){
974  NNEW(hcel,double,*nef);
975  FORTRAN(calchcel,(vel,cocon,ncocon,ielmatf,ntmat_,mi,
976  hcel,nef));
977  }
978 
979  /* calculating the stress and/or the heat flux at the
980  element centers */
981 
982  if((isti==1)||(iqfx==1)){
983  FORTRAN(calcstressheatflux,(sti,umel,gradvel,qfx,hcel,
984  gradtel,nef,&isti,&iqfx,mi));
985  if(iqfx==1)SFREE(hcel);
986  }
987 
988  /* extrapolating the stresses */
989 
990  if((strcmp1(&filab[3306],"SF ")==0)||
991  (strcmp1(&filab[3480],"SVF ")==0)){
992  nfield=6;
993  ndim=6;
994  if((*norien>0)&&
995  ((strcmp1(&filab[3311],"L")==0)||(strcmp1(&filab[3485],"L")==0))){
996  iorienglob=1;
997  }else{
998  iorienglob=0;
999  }
1000  strcpy1(&cflag[0],&filab[2962],1);
1001  NNEW(stn,double,6**nk);
1002  FORTRAN(extrapolate,(sti,stn,ipkonf,inum,kon,lakonf,
1003  &nfield,nk,nef,mi,&ndim,orab,ielorienf,co,&iorienglob,
1004  cflag,nelemload,nload,nodeboun,nboun,ndirboun,
1005  vold,ithermal,&force,&icfdout,ielmatf,thicke,filab));
1006  }
1007 
1008  /* extrapolating the heat flow */
1009 
1010 
1011  if(strcmp1(&filab[3393],"HFLF")==0){
1012  nfield=3;
1013  ndim=3;
1014  if((*norien>0)&&(strcmp1(&filab[3398],"L")==0)){
1015  iorienglob=1;
1016  }else{
1017  iorienglob=0;
1018  }
1019  strcpy1(&cflag[0],&filab[3049],1);
1020  NNEW(qfn,double,3**nk);
1021  FORTRAN(extrapolate,(qfx,qfn,ipkonf,inum,kon,lakonf,
1022  &nfield,nk,nef,mi,&ndim,orab,ielorienf,co,&iorienglob,
1023  cflag,nelemload,nload,nodeboun,nboun,ndirboun,
1024  vold,ithermal,&force,&icfdout,ielmatf,thicke,filab));
1025  }
1026 
1027  /* extrapolating the facial values of the static temperature
1028  and/or the velocity and/or the static pressure to the nodes */
1029 
1030  FORTRAN(extrapolatefluid,(nk,iponofa,inofa,inum,vfa,vold,ielfa,
1031  ithermal));
1032 
1033  /* storing the results in dat-format */
1034 
1035  ptimef=ttimef+timef;
1036  FORTRAN(printoutfluid,(set,nset,istartset,iendset,ialset,nprint,
1037  prlab,prset,v,t1,fn,ipkonf,lakonf,sti,eei,
1038  xstate,ener,mi,&nstate_,ithermal,co,kon,qfx,
1039  &ptimef,trab,inotr,ntrans,orab,ielorienf,
1040  norien,nk,nef,inum,filab,vold,ielmatf,
1041  thicke,eme,vcontu,physcon,nactdoh,
1042  ielpropf,prop));
1043 
1044  /* storing the results in frd-format */
1045 
1046  FORTRAN(frdfluid,(co,nk,kon,ipkonf,lakonf,nef,vold,&kode,&timef,ielmatf,
1047  matname,filab,inum,ntrans,inotr,trab,mi,istep,
1048  stn,qfn,nactdohinv));
1049 
1050  if((strcmp1(&filab[3306],"SF ")==0)||
1051  (strcmp1(&filab[3480],"SVF ")==0)){SFREE(stn);}
1052  if(strcmp1(&filab[3393],"HFLF")==0){SFREE(qfn);}
1053 
1054  }
1055 
1056  if(iit==jmax[1]){
1057  printf("*INFO: maximum number of fluid increments reached\n\n");
1058  FORTRAN(stop,());
1059  }
1060  if(last==1){FORTRAN(stop,());}
1061 
1062  memcpy(&veloo[0],&velo[0],sizeof(double)*6**nef);
1063  memcpy(&velo[0],&vel[0],sizeof(double)*6**nef);
1064 
1065  }while(1);
1066 
1067  FORTRAN(closefilefluid,());
1068 
1069  SFREE(irow);SFREE(ia);SFREE(icol);SFREE(jq);SFREE(ja);
1070 //SFREE(nactdoh);
1071 
1072  SFREE(coel);SFREE(cosa);SFREE(xxn);SFREE(xxi);SFREE(xle);SFREE(xlen);
1073  SFREE(xlet);SFREE(cofa);SFREE(area);SFREE(xrlfa);SFREE(volume);
1074  SFREE(cosb);SFREE(xxj);
1075 
1076  SFREE(ifabou);SFREE(umfa);SFREE(umel);
1077 
1078  SFREE(gradvel);SFREE(gradvfa);SFREE(au);SFREE(ad);SFREE(b);SFREE(advfa);
1079  SFREE(ap);SFREE(bp);SFREE(gradpel);SFREE(rwork);SFREE(aua);
1080  SFREE(hfa);SFREE(hel);SFREE(auv);SFREE(adv);SFREE(bv);SFREE(sel);
1081 
1082  if(*ithermal>0){
1083  SFREE(gradtel);SFREE(gradtfa);SFREE(hcfa);SFREE(cpel);SFREE(cpfa);
1084  }
1085 
1086  SFREE(inum);SFREE(v);SFREE(velo);SFREE(veloo);
1087 
1088  SFREE(iponofa);SFREE(inofa);
1089 
1090  if(*nbody>0) SFREE(body);
1091 
1092  *ithermal=ithermalref;
1093 
1094  return;
1095 
1096 }
#define ITGFORMAT
Definition: CalculiX.h:52
subroutine calchcel(vel, cocon, ncocon, ielmat, ntmat_, mi, hcel, nef)
Definition: calchcel.f:19
subroutine extrapolate_gradvel(nface, ielfa, xrlfa, gradvel, gradvfa)
Definition: extrapolate_gradvel.f:19
subroutine calcgradvel(ne, lakon, ipnei, vfa, area, xxn, gradvel, neifa, volume)
Definition: calcgradvel.f:19
subroutine extrapolate_pel(nface, ielfa, xrlfa, vel, vfa, ifabou, xboun, nef)
Definition: extrapolate_pel.f:19
subroutine calcgradpel(nef, lakon, ipnei, vfa, area, xxn, gradpel, neifa, volume)
Definition: calcgradpel.f:19
subroutine calcbody(nef, body, ipobody, ibody, xbody, coel, vel, lakon, nactdohinv)
Definition: calcbody.f:19
subroutine correctvfa(nface, ielfa, area, vfa, ap, bp, xxn, ifabou, ipnei, nef, neifa, hfa, vel, xboun)
Definition: correctvfa.f:19
subroutine initialcfd(nef, ipkonf, kon, lakonf, co, coel, cofa, nface, ielfa, area, ipnei, neiel, xxn, xxi, xle, xlen, xlet, xrlfa, cosa, volume, neifa, xxj, cosb, vel, dmin)
Definition: initialcfd.f:19
subroutine matvec(n, x, y, nelt, ia, ja, a, isym)
Definition: matvec.f:26
subroutine inicalcbody(nef, body, ipobody, ibody, xbody, coel, vel, lakon, nactdohinv, icent)
Definition: inicalcbody.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))
subroutine extrapolate_gradtel(nface, ielfa, xrlfa, gradtel, gradtfa)
Definition: extrapolate_gradtel.f:19
subroutine calccpel(nef, vel, shcon, nshcon, ielmat, ntmat_, mi, cpel, physcon)
Definition: calccpel.f:20
subroutine calcumel(nef, vel, shcon, nshcon, ielmat, ntmat_, ithermal, mi, umel)
Definition: calcumel.f:19
subroutine calcrhofa(nface, vfa, rhcon, nrhcon, ielmat, ntmat_, ithermal, mi, ielfa)
Definition: calcrhofa.f:19
subroutine extrapolatefluid(nk, iponofa, inofa, inum, vfa, v, ielfa, ithermal)
Definition: extrapolatefluid.f:19
ITG strcpy1(char *s1, const char *s2, ITG length)
Definition: strcpy1.c:24
subroutine printoutfluid(set, nset, istartset, iendset, ialset, nprint, prlab, prset, v, t1, fn, ipkonf, lakonf, sti, eei, xstate, ener, mi, nstate_, ithermal, co, kon, qfx, ttime, trab, inotr, ntrans, orab, ielorienf, norien, nk, ne, inum, filab, vold, ielmatf, thicke, eme, vcontu, physcon, nactdoh, ielpropf, prop)
Definition: printoutfluid.f:19
ITG strcmp1(const char *s1, const char *s2)
Definition: strcmp1.c:24
subroutine mafillp(nef, lakon, ipnei, neifa, neiel, vfa, area, advfa, xlet, cosa, volume, au, ad, jq, irow, ap, ielfa, ifabou, xle, b, xxn, compressible, neq, nzs, hfa, gradpel, bp, xxi, neij, xlen, cosb)
Definition: mafillp.f:19
void mastructf(ITG *nk, ITG *kon, ITG *ipkon, char *lakon, ITG *ne, ITG *icol, ITG *jq, ITG **mast1p, ITG **irowp, ITG *isolver, ITG *neq, ITG *ipointer, ITG *nzs, ITG *ipnei, ITG *ineiel, ITG *mi)
Definition: mastructf.c:27
subroutine openfilefluid(jobname)
Definition: openfilefluid.f:19
ITG getSystemCPUs()
Definition: getSystemCPUs.c:40
#define DMEMSET(a, b, c, d)
Definition: CalculiX.h:45
subroutine stop()
Definition: stop.f:19
subroutine dslugm(N, B, X, NELT, IA, JA, A, ISYM, NSAVE, ITOL, TOL, ITMAX, ITER, ERR, IERR, IUNIT, RWORK, LENW, IWORK, LENIW)
Definition: dslugm.f:2
subroutine rhsp(nef, lakon, nactdoh, ipnei, neifa, neiel, vfa, area, advfa, xlet, cosa, volume, au, ad, jq, irow, ap, ielfa, ifabou, xle, b, xxn, compressible, neq, nzs, hfa, bp, neij, xxi, gradpel, xlen)
Definition: rhsp.f:19
subroutine convert2slapcol(au, ad, irow, ia, jq, ja, nzs, neq, aua)
Definition: convert2slapcol.f:25
subroutine applympc(nface, ielfa, is, ie, ifabou, ipompc, vfa, coefmpc, nodempc, ipnei, neifa, labmpc, xbounact, nactdoh)
Definition: applympc.f:19
subroutine extrapolate(yi, yn, ipkon, inum, kon, lakon, nfield, nk, ne, mi, ndim, orab, ielorien, co, iorienloc, cflag, vold, force, ielmat, thicke, ielprop, prop)
Definition: extrapolate.f:19
subroutine closefilefluid()
Definition: closefilefluid.f:19
subroutine calchcfa(nface, vfa, cocon, ncocon, ielmat, ntmat_, mi, ielfa, hcfa)
Definition: calchcfa.f:19
subroutine calcstressheatflux(sti, umel, gradvel, qfx, hcel, gradtel, nef, isti, iqfx, mi)
Definition: calcstressheatflux.f:19
subroutine calcumfa(nface, vfa, shcon, nshcon, ielmat, ntmat_, ithermal, mi, ielfa, umfa)
Definition: calcumfa.f:19
subroutine predgmres(n, b, x, nelt, ia, ja, a, isym, itol, tol, itmax, iter, err, ierr, iunit, sb, sx, rgwk, lrgw, igwk, ligw, rwork, iwork)
Definition: predgmres.f:26
subroutine frdfluid(co, nk, kon, ipkonf, lakonf, nef, vold, kode, time, ielmat, matname, filab, inum, ntrans, inotr, trab, mi, istep, stn, qfn, nactdohinv)
Definition: frdfluid.f:19
#define RENEW(a, b, c)
Definition: CalculiX.h:40
#define SFREE(a)
Definition: CalculiX.h:41
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 complete_hel(neq, bv, hel, adv, auv, jq, irow, nzs)
Definition: complete_hel.f:25
subroutine mafillt(nef, ipnei, neifa, neiel, vfa, xxn, area, au, ad, jq, irow, nzs, b, vel, umel, xlet, xle, gradtfa, xxi, body, volume, compressible, ielfa, lakon, ifabou, nbody, neq, dtimef, velo, veloo, cpfa, hcfa, cpel, gradvel, xload, gammat, xrlfa, xxj, nactdohinv, a1, a2, a3)
Definition: mafillt.f:19
static ITG num_cpus
Definition: compfluid.c:37
subroutine correctvel(hel, adv, vfa, ipnei, area, bv, xxn, neifa, lakon, ne, neq)
Definition: correctvel.f:19
subroutine calcrhoel(nef, vel, rhcon, nrhcon, ielmat, ntmat_, ithermal, mi)
Definition: calcrhoel.f:19
subroutine extrapolate_vel(nface, ielfa, xrlfa, vel, vfa, ifabou, xboun, ipnei, nef)
Definition: extrapolate_vel.f:19
subroutine applyboun(ifaext, nfaext, ielfa, ikboun, ilboun, nboun, typeboun, nelemload, nload, sideload, isolidsurf, nsolidsurf, ifabou, nfabou, nface, nodeboun, ndirboun, ikmpc, ilmpc, labmpc, nmpc, nactdohinv)
Definition: applyboun.f:19
subroutine extrapolate_tel(nface, ielfa, xrlfa, vel, vfa, ifabou, xboun, ipnei, nef)
Definition: extrapolate_tel.f:19
subroutine calcguesstincf(nface, dmin, vfa, umfa, tincfguess)
Definition: calcguesstincf.f:19
subroutine cataloguenodes(iponofa, inofa, ifreefa, ielfa, ifabou, ipkon, kon, lakon, nface, nk)
Definition: cataloguenodes.f:19
#define ITG
Definition: CalculiX.h:51
subroutine msolve(n, r, z, nelt, ia, ja, a, isym, rwork, iwork)
Definition: msolve.f:21
subroutine calcgradtel(nef, lakon, ipnei, vfa, area, xxn, gradtel, neifa, volume)
Definition: calcgradtel.f:19
subroutine calccpfa(nface, vfa, shcon, nshcon, ielmat, ntmat_, mi, ielfa, cpfa, physcon)
Definition: calccpfa.f:19
subroutine extrapolate_ad_h(nface, ielfa, xrlfa, adv, advfa, hel, hfa)
Definition: extrapolate_ad_h.f:19
#define NNEW(a, b, c)
Definition: CalculiX.h:39
subroutine mafillv(nef, ipnei, neifa, neiel, vfa, xxn, area, auv, adv, jq, irow, nzs, bv, vel, cosa, umfa, xlet, xle, gradvfa, xxi, body, volume, compressible, ielfa, lakon, ifabou, nbody, neq, dtimef, velo, veloo, sel, xrlfa, gamma, xxj, nactdohinv, a1, a2, a3)
Definition: mafillv.f:19

Variable Documentation

ITG num_cpus
static

Definition at line 37 of file compfluid.c.

Hosted by OpenAircraft.com, (Michigan UAV, LLC)