/* CalculiX - A 3-dimensional finite element program */ /* Copyright (C) 1998-2007 Guido Dhondt */ /* This program is free software; you can redistribute it and/or */ /* modify it under the terms of the GNU General Public License as */ /* published by the Free Software Foundation(version 2); */ /* */ /* This program is distributed in the hope that it will be useful, */ /* but WITHOUT ANY WARRANTY; without even the implied warranty of */ /* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the */ /* GNU General Public License for more details. */ /* You should have received a copy of the GNU General Public License */ /* along with this program; if not, write to the Free Software */ /* Foundation, Inc., 675 Mass Ave, Cambridge, MA 02139, USA. */ /* calculating the initial acceleration at the start of the step for dynamic calculations */ if((*nmethod==4)&&(*ithermal!=2)){ bet=(1.-*alpha)*(1.-*alpha)/4.; gam=0.5-*alpha; /* calculating the stiffness and mass matrix the stress must be determined to calculate the stiffness matrix*/ reltime=0.; time=0.; dtime=0.; FORTRAN(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)); time=0.; dtime=1.; /* updating the nonlinear mpc's (also affects the boundary conditions through the nonhomogeneous part of the mpc's) if contact arises the number of MPC's can also change */ cam[0]=0.;cam[1]=0.; if(*ithermal>1){radflowload(itg,ieg,&ntg,&ntr,&ntm, ac,bc,nload,sideload,nelemload,xloadact,lakon,ipiv,ntmat_,vold, shcon,nshcon,ipkon,kon,co,pmid,e1,e2,e3,iptri, kontri,&ntri,nloadtr,tarea,tenv,physcon,erad,fij, dist,idist,area,nflow,ikboun,xboun,nboun,ithermal,&iinc,&iit, cs,mcs,inocs,&ntrit,nk,fenv,istep,&dtime,ttime,&time,ilboun, ikforc,ilforc,xforcact,nforc,cam,ielmat,&nteq,prop,ielprop, nactdog,nacteq,nodeboun,ndirboun,voldgas,&network,rhcon, nrhcon,ipobody,ibody,xbodyact,nbody,iviewfile,jobnamef,ctrl);} if(ncont!=0){ *ne=ne0;*nkon=nkon0; contact(&ncont,ntie,tieset,nset,set,istartset,iendset, ialset,itietri,lakon,ipkon,kon,koncont,ne,cg,straight,nkon, co,vold,ielmat,cs); icascade=1; } if(icascade==2){ *nmpc=nmpcref;memmpc_=memmpcref_;mpcfree=mpcfreeref; RENEW(nodempc,int,3*memmpcref_); for(k=0;k<3*memmpcref_;k++){nodempc[k]=nodempcref[k];} RENEW(coefmpc,double,memmpcref_); for(k=0;k0) remastruct(ipompc,&coefmpc,&nodempc,nmpc, &mpcfree,nodeboun,ndirboun,nboun,ikmpc,ilmpc,ikboun,ilboun, labmpc,nk,&memmpc_,&icascade,&maxlenmpc, kon,ipkon,lakon,ne,nnn,nactdof,icol,jq,&irow,isolver, neq,nzs,nmethod,&f,&fext,&b,&aux2,&fini,&fextini, &adb,&aub,ithermal,iperturb); iout=-1; ielas=1; fn=NNEW(double,4**nk); stx=NNEW(double,6**mint_**ne); if(*ithermal>1) qfx=NNEW(double,3**mint_**ne); FORTRAN(results,(co,nk,kon,ipkon,lakon,ne,vold,stn,inum,stx, elcon,nelcon,rhcon,nrhcon,alcon,nalcon,alzero,ielmat, ielorien,norien,orab,ntmat_,t0,t1old,ithermal, prestr,iprestr,filab,eme,een,iperturb, f,fn,nactdof,&iout,qa,vold,b,nodeboun, ndirboun,xbounold,nboun,ipompc, nodempc,coefmpc,labmpc,nmpc,nmethod,cam,&neq[1],veold,accold,&bet, &gam,&dtime,&time,ttime,plicon,nplicon,plkcon,nplkcon, xstateini,xstiff,xstate,npmat_,epn,matname,mint_,&ielas,&icmd, ncmat_,nstate_,sti,vini,ikboun,ilboun,ener,enern,sti,xstaten, eei,enerini,cocon,ncocon,set,nset,istartset,iendset, ialset,nprint,prlab,prset,qfx,qfn,trab,inotr,ntrans,fmpc, nelemload,nload)); free(fn);free(stx);if(*ithermal>1)free(qfx); iout=0; ielas=0; reltime=0.; time=0.; dtime=0.; if(*iexpl==0){intscheme=1;} /* in mafillsm the stiffness and mass matrix are computed; The primary aim is to calculate the mass matrix (not lumped for an implicit dynamic calculation, lumped for an explicit dynamic calculation). However: - for an implicit calculation the mass matrix is "doped" with a small amount of stiffness matrix, therefore the calculation of the stiffness matrix is needed. - for an explicit calculation the stiffness matrix is not needed at all. Since the calculation of the mass matrix alone is not possible in mafillsm, the determination of the stiffness matrix is taken as unavoidable "ballast". */ ad=NNEW(double,neq[1]); au=NNEW(double,nzs[1]); FORTRAN(mafillsm,(co,nk,kon,ipkon,lakon,ne,nodeboun,ndirboun,xbounold,nboun, ipompc,nodempc,coefmpc,nmpc,nodeforc,ndirforc,xforcact, nforc,nelemload,sideload,xloadact,nload,xbodyact,ipobody, nbody,cgr,ad,au,fext,nactdof,icol,jq,irow,neq,nzl, &nmethodact,ikmpc,ilmpc,ikboun,ilboun, elcon,nelcon,rhcon,nrhcon,alcon,nalcon,alzero, ielmat,ielorien,norien,orab,ntmat_, t0,t1act,ithermal,prestr,iprestr,vold,iperturb,sti, nzs,stx,adb,aub,iexpl,plicon,nplicon,plkcon,nplkcon, xstiff,npmat_,&dtime,matname,mint_, ncmat_,mass,&stiffness,&buckling,&rhsi,&intscheme, physcon,shcon,nshcon,cocon,ncocon,ttime,&time,istep,&iinc, &coriolis,ibody)); if(nmethodact==0){ /* error occurred in mafill: storing the geometry in frd format */ ++*kode; ipneigh=NNEW(int,*nk);neigh=NNEW(int,40**ne); FORTRAN(out,(co,nk,kon,ipkon,lakon,ne,v,stn,inum,nmethod,kode,filab, een,t1,fn,ttime,epn,ielmat,matname,enern,xstaten,nstate_,istep, &iinc,iperturb,ener,mint_,output,ithermal,qfn,&mode,&noddiam, trab,inotr,ntrans,orab,ielorien,norien,description, ipneigh,neigh,sti));free(ipneigh);free(neigh); FORTRAN(stop,()); } /* calculating the acceleration at the start of the step. This can be different from the acceleration at the end of the last step due to a discontinuous loading increase */ /* reltime=0.; time=0.; dtime=0.; FORTRAN(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));*/ /* determining the external loading vector */ /* FORTRAN(rhs,(co,nk,kon,ipkon,lakon,ne, ipompc,nodempc,coefmpc,nmpc,nodeforc,ndirforc,xforcact, nforc,nelemload,sideload,xloadact,nload,xbodyact,ipobody, nbody,cgr,fext,nactdof,&neq[1], &nmethodact,ikmpc,ilmpc, elcon,nelcon,rhcon,nrhcon,alcon,nalcon,alzero, ielmat,ielorien,norien,orab,ntmat_, t0,t1act,ithermal,iprestr,vold,iperturb, iexpl,plicon,nplicon,plkcon,nplkcon, npmat_,ttime,&time,istep,&iinc,&dtime,physcon,ibody));*/ /* mass x acceleration = f(external)-f(internal) only for the mechanical loading*/ for(k=0;k