1! 2! CalculiX - A 3-dimensional finite element program 3! Copyright (C) 1998-2021 Guido Dhondt 4! 5! This program is free software; you can redistribute it and/or 6! modify it under the terms of the GNU General Public License as 7! published by the Free Software Foundation(version 2); 8! 9! 10! This program is distributed in the hope that it will be useful, 11! but WITHOUT ANY WARRANTY; without even the implied warranty of 12! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the 13! GNU General Public License for more details. 14! 15! You should have received a copy of the GNU General Public License 16! along with this program; if not, write to the Free Software 17! Foundation, Inc., 675 Mass Ave, Cambridge, MA 02139, USA. 18! 19 subroutine restartread(istep,nset,nload,nforc, nboun,nk,ne, 20 & nmpc,nalset,nmat,ntmat_,npmat_,norien,nam,nprint,mi, 21 & ntrans,ncs_,namtot,ncmat_,mpcfree,maxlenmpc, 22 & ne1d,ne2d,nflow,nlabel,iplas, 23 & nkon,ithermal,nmethod,iperturb,nstate_,nener,set,istartset, 24 & iendset,ialset,co,kon,ipkon,lakon,nodeboun,ndirboun,iamboun, 25 & xboun,ikboun,ilboun,ipompc,nodempc,coefmpc,labmpc,ikmpc,ilmpc, 26 & nodeforc,ndirforc,iamforc,xforc,ikforc,ilforc,nelemload,iamload, 27 & sideload,xload,elcon,nelcon,rhcon,nrhcon, 28 & alcon,nalcon,alzero,plicon,nplicon,plkcon,nplkcon,orname,orab, 29 & ielorien,trab,inotr,amname,amta,namta,t0,t1,iamt1,veold, 30 & ielmat,matname,prlab,prset,filab,vold,nodebounold, 31 & ndirbounold,xbounold,xforcold,xloadold,t1old,eme, 32 & iponor,xnor,knor,thicke,offset,iponoel,inoel,rig, 33 & shcon,nshcon,cocon,ncocon,ics,sti, 34 & ener,xstate,jobnamec,infree,irestartstep,prestr,iprestr, 35 & cbody,ibody,xbody,nbody,xbodyold,ttime,qaold,cs,mcs, 36 & output,physcon,ctrl,typeboun,fmpc,tieset,ntie,tietol,nslavs, 37 & t0g,t1g,nprop,ielprop,prop,mortar,nintpoint,ifacecount, 38 & islavsurf,pslavsurf,clearini,irstrt,vel,nef,velo,veloo, 39 & ne2boun,heading,network) 40! 41 implicit none 42! 43 character*1 typeboun(*) 44 character*4 output 45 character*6 prlab(*) 46 character*8 lakon(*) 47 character*20 labmpc(*),sideload(*) 48 character*66 heading(*) 49 character*80 orname(*),amname(*),matname(*),version 50 character*81 set(*),prset(*),tieset(*),cbody(*) 51 character*87 filab(*) 52 character*132 fnrstrt,jobnamec(*) 53! 54 integer istep,nset,nload,nforc,nboun,nk,ne,nmpc,nalset,nmat, 55 & ntmat_,npmat_,norien,nam,nprint,mi(*),ntrans,ncs_,nheading_, 56 & namtot,ncmat_,mpcfree,ne1d,ne2d,nflow,nlabel,iplas,nkon, 57 & ithermal(*),nmethod,iperturb(*),nstate_,istartset(*),iendset(*), 58 & ialset(*),kon(*),ipkon(*),nodeboun(*),ndirboun(*),iamboun(*), 59 & ikboun(*),ilboun(*),ipompc(*),nodempc(*),ikmpc(*),ilmpc(*), 60 & nodeforc(*),ndirforc(*),iamforc(*),ikforc(*),ilforc(*), 61 & nelemload(*),iamload(*),nelcon(*),mt,nprop,ielprop(*), 62 & nrhcon(*),nalcon(*),nplicon(*),nplkcon(*),ielorien(*), 63 & inotr(*),mortar,nintpoint,ifacecount,islavsurf(*), 64 & namta(*),iamt1(*),ielmat(*),nodebounold(*),ndirbounold(*), 65 & iponor(*),knor(*),iponoel(*),inoel(*),rig(*), 66 & nshcon(*),ncocon(*),ics(*),infree(*),i,ipos, 67 & nener,irestartstep,istat,iprestr,irstrt(*), 68 & maxlenmpc,mcs,mpcend,ntie,ibody(*),nbody,nslavs,nef, 69 & ne2boun(*),memmpc_,network 70! 71 real*8 co(*),xboun(*),coefmpc(*),xforc(*),xload(*),elcon(*), 72 & rhcon(*),alcon(*),alzero(*),plicon(*),plkcon(*),orab(*), 73 & trab(*),amta(*),t0(*),t1(*),veold(*),tietol(*),pslavsurf(*), 74 & vold(*),xbounold(*),xforcold(*),xloadold(*),t1old(*),eme(*), 75 & xnor(*),thicke(*),offset(*),t0g(*),t1g(*),clearini(*), 76 & shcon(*),cocon(*),sti(*),ener(*),xstate(*),prestr(*),ttime, 77 & qaold(2),physcon(*),ctrl(*),cs(*),fmpc(*),xbody(*), 78 & xbodyold(*),prop(*),vel(*),velo(*),veloo(*) 79! 80 ipos=index(jobnamec(1),char(0)) 81 fnrstrt(1:ipos-1)=jobnamec(1)(1:ipos-1) 82 fnrstrt(ipos:ipos+3)=".rin" 83 do i=ipos+4,132 84 fnrstrt(i:i)=' ' 85 enddo 86! 87 open(15,file=fnrstrt,ACCESS='SEQUENTIAL',FORM='UNFORMATTED', 88 & err=15) 89! 90 do 91! 92 read(15,iostat=istat) version 93! 94 if(istat.lt.0) then 95 if(irestartstep.eq.0) then 96! 97! reading the last step 98! 99 irestartstep=istep 100 close(15) 101 open(15,file=fnrstrt,ACCESS='SEQUENTIAL', 102 & FORM='UNFORMATTED',err=15) 103 read(15) version 104 else 105 write(*,*) '*ERROR reading *RESTART,READ: requested step' 106 write(*,*) ' is not in the restart file' 107 call exit(201) 108 endif 109 endif 110 read(15)istep 111! 112! set size 113! 114 read(15)nset 115 read(15)nalset 116! 117! load size 118! 119 read(15)nload 120 read(15)nbody 121 read(15)nforc 122 read(15)nboun 123 read(15)nflow 124! 125! mesh size 126! 127 read(15)nk 128 read(15)ne 129 read(15)nef 130 read(15)nkon 131 read(15)(mi(i),i=1,3) 132 mt=mi(2)+1 133! 134! constraint size 135! 136 read(15)nmpc 137 read(15)mpcend 138 read(15)maxlenmpc 139 read(15)memmpc_ 140! 141! material size 142! 143 read(15)nmat 144 read(15)ntmat_ 145 read(15)npmat_ 146 read(15)ncmat_ 147! 148! property info 149! 150 read(15)nprop 151! 152! transformation size 153! 154 read(15)norien 155 read(15)ntrans 156! 157! amplitude size 158! 159 read(15)nam 160 read(15)namtot 161! 162! print size 163! 164 read(15)nprint 165 read(15)nlabel 166! 167! tie size 168! 169 read(15)ntie 170! 171! cyclic symmetry size 172! 173 read(15)ncs_ 174 read(15)mcs 175! 176! 1d and 2d element size 177! 178 read(15)ne1d 179 read(15)ne2d 180 read(15)(infree(i),i=1,4) 181! 182! procedure info 183! 184 read(15)nmethod 185 read(15)(iperturb(i),i=1,2) 186 read(15)nener 187 read(15)iplas 188 read(15)ithermal(1) 189 read(15)nstate_ 190 read(15)nslavs 191 read(15)iprestr 192 read(15)mortar 193 if(mortar.eq.1) then 194 read(15)ifacecount 195 read(15)nintpoint 196 endif 197 read(15)nheading_ 198 read(15)network 199! 200 if(istep.eq.irestartstep) exit 201! 202! skipping the next entries until the requested step is found 203! 204 call skip(nset,nalset,nload,nbody,nforc,nboun,nk,ne,nkon, 205 & mi(1),nmpc,mpcend,nmat,ntmat_,npmat_,ncmat_,norien,ntrans, 206 & nam,nprint,nlabel,ncs_,ne1d,ne2d,infree,nmethod, 207 & iperturb,nener,ithermal,nstate_,iprestr,mcs,ntie, 208 & nslavs,nprop,mortar,ifacecount,nintpoint,nef,nheading_) 209! 210 enddo 211! 212! sets 213! 214 read(15)(set(i),i=1,nset) 215 read(15)(istartset(i),i=1,nset) 216 read(15)(iendset(i),i=1,nset) 217 do i=1,nalset 218 read(15)ialset(i) 219 enddo 220! 221! header lines 222! 223 read(15)(heading(i),i=1,nheading_) 224! 225! mesh 226! 227 read(15)(co(i),i=1,3*nk) 228 read(15)(kon(i),i=1,nkon) 229 read(15)(ipkon(i),i=1,ne) 230 read(15)(lakon(i),i=1,ne) 231! 232! single point constraints 233! 234 read(15)(nodeboun(i),i=1,nboun) 235 read(15)(ndirboun(i),i=1,nboun) 236 read(15)(typeboun(i),i=1,nboun) 237 read(15)(xboun(i),i=1,nboun) 238 read(15)(ikboun(i),i=1,nboun) 239 read(15)(ilboun(i),i=1,nboun) 240 if(nam.gt.0) read(15)(iamboun(i),i=1,nboun) 241 read(15)(nodebounold(i),i=1,nboun) 242 read(15)(ndirbounold(i),i=1,nboun) 243 read(15)(xbounold(i),i=1,nboun) 244! 245! multiple point constraints 246! 247 read(15)(ipompc(i),i=1,nmpc) 248 read(15)(labmpc(i),i=1,nmpc) 249 read(15)(ikmpc(i),i=1,nmpc) 250 read(15)(ilmpc(i),i=1,nmpc) 251 read(15)(fmpc(i),i=1,nmpc) 252 read(15)(nodempc(i),i=1,3*mpcend) 253 read(15)(coefmpc(i),i=1,mpcend) 254 mpcfree=mpcend+1 255! 256! point forces 257! 258 read(15)(nodeforc(i),i=1,2*nforc) 259 read(15)(ndirforc(i),i=1,nforc) 260 read(15)(xforc(i),i=1,nforc) 261 read(15)(ikforc(i),i=1,nforc) 262 read(15)(ilforc(i),i=1,nforc) 263 if(nam.gt.0) read(15)(iamforc(i),i=1,nforc) 264 read(15)(xforcold(i),i=1,nforc) 265! 266! distributed loads 267! 268 read(15)(nelemload(i),i=1,2*nload) 269 read(15)(sideload(i),i=1,nload) 270 read(15)(xload(i),i=1,2*nload) 271 if(nam.gt.0) read(15)(iamload(i),i=1,2*nload) 272 read(15)(xloadold(i),i=1,2*nload) 273 read(15)(cbody(i),i=1,nbody) 274 read(15)(ibody(i),i=1,3*nbody) 275 read(15)(xbody(i),i=1,7*nbody) 276 read(15)(xbodyold(i),i=1,7*nbody) 277! 278! prestress 279! 280 if(iprestr.gt.0) read(15) (prestr(i),i=1,6*mi(1)*ne) 281! 282! labels 283! 284 read(15)(prlab(i),i=1,nprint) 285 read(15)(prset(i),i=1,nprint) 286 read(15)(filab(i),i=1,nlabel) 287! 288! elastic constants 289! 290 read(15)(elcon(i),i=1,(ncmat_+1)*ntmat_*nmat) 291 read(15)(nelcon(i),i=1,2*nmat) 292! 293! density 294! 295 read(15)(rhcon(i),i=1,2*ntmat_*nmat) 296 read(15)(nrhcon(i),i=1,nmat) 297! 298! specific heat 299! 300 read(15)(shcon(i),i=1,4*ntmat_*nmat) 301 read(15)(nshcon(i),i=1,nmat) 302! 303! conductivity 304! 305 read(15)(cocon(i),i=1,7*ntmat_*nmat) 306 read(15)(ncocon(i),i=1,2*nmat) 307! 308! expansion coefficients 309! 310 read(15)(alcon(i),i=1,7*ntmat_*nmat) 311 read(15)(nalcon(i),i=1,2*nmat) 312 read(15)(alzero(i),i=1,nmat) 313! 314! physical constants 315! 316 read(15)(physcon(i),i=1,10) 317! 318! plastic data 319! 320 if(npmat_.ne.0)then 321 read(15)(plicon(i),i=1,(2*npmat_+1)*ntmat_*nmat) 322 read(15)(nplicon(i),i=1,(ntmat_+1)*nmat) 323 read(15)(plkcon(i),i=1,(2*npmat_+1)*ntmat_*nmat) 324 read(15)(nplkcon(i),i=1,(ntmat_+1)*nmat) 325 endif 326! 327! material orientation 328! 329 if(norien.ne.0)then 330 read(15)(orname(i),i=1,norien) 331 read(15)(orab(i),i=1,7*norien) 332 read(15)(ielorien(i),i=1,mi(3)*ne) 333 endif 334! 335! fluid section properties 336! 337 if(nprop.ne.0) then 338 read(15)(ielprop(i),i=1,ne) 339 read(15)(prop(i),i=1,nprop) 340 endif 341! 342! transformations 343! 344 if(ntrans.ne.0)then 345 read(15)(trab(i),i=1,7*ntrans) 346 read(15)(inotr(i),i=1,2*nk) 347 endif 348! 349! amplitudes 350! 351 if(nam.gt.0)then 352 read(15)(amname(i),i=1,nam) 353 read(15)(namta(i),i=1,3*nam-1) 354 read(15) namta(3*nam) 355 read(15)(amta(i),i=1,2*namta(3*nam-1)) 356 endif 357! 358! temperatures 359! 360 if(ithermal(1).gt.0)then 361 read(15)(t0(i),i=1,nk) 362 read(15)(t1(i),i=1,nk) 363 if((ne1d.gt.0).or.(ne2d.gt.0))then 364 read(15)(t0g(i),i=1,2*nk) 365 read(15)(t1g(i),i=1,2*nk) 366 endif 367 if(nam.gt.0) read(15)(iamt1(i),i=1,nk) 368 read(15)(t1old(i),i=1,nk) 369 endif 370! 371! materials 372! 373 read(15)(matname(i),i=1,nmat) 374 read(15)(ielmat(i),i=1,mi(3)*ne) 375! 376! temperature, displacement, static pressure, velocity and acceleration 377! 378 read(15)(vold(i),i=1,mt*nk) 379 if((nmethod.eq.4).or.((nmethod.eq.1).and.(iperturb(1).ge.2))) then 380 read(15)(veold(i),i=1,mt*nk) 381 endif 382! 383! CFD results at the element centers 384! 385 if(nef.gt.0) then 386 read(15)(vel(i),i=1,8*nef) 387 read(15)(velo(i),i=1,8*nef) 388 read(15)(veloo(i),i=1,8*nef) 389 endif 390! 391! 1d and 2d elements 392! 393 if((ne1d.gt.0).or.(ne2d.gt.0))then 394 read(15)(iponor(i),i=1,2*nkon) 395 read(15)(xnor(i),i=1,infree(1)) 396 read(15)(knor(i),i=1,infree(2)) 397 read(15)(thicke(i),i=1,mi(3)*nkon) 398 read(15)(offset(i),i=1,2*ne) 399 read(15)(iponoel(i),i=1,infree(4)) 400 read(15)(inoel(i),i=1,3*(infree(3)-1)) 401 read(15)(rig(i),i=1,infree(4)) 402 read(15)(ne2boun(i),i=1,2*infree(4)) 403 endif 404! 405! tie constraints 406! 407 if(ntie.gt.0) then 408 read(15)(tieset(i),i=1,3*ntie) 409 read(15)(tietol(i),i=1,3*ntie) 410 endif 411! 412! cyclic symmetry 413! 414 if(ncs_.gt.0)then 415 read(15)(ics(i),i=1,ncs_) 416 endif 417 if(mcs.gt.0) then 418 read(15)(cs(i),i=1,17*mcs) 419 endif 420! 421! integration point variables 422! 423 read(15)(sti(i),i=1,6*mi(1)*ne) 424 read(15)(eme(i),i=1,6*mi(1)*ne) 425 if(nener.eq.1) then 426 read(15)(ener(i),i=1,mi(1)*ne) 427 endif 428 if(nstate_.gt.0)then 429 if(mortar.eq.0) then 430 read(15)(xstate(i),i=1,nstate_*mi(1)*(ne+nslavs)) 431 elseif(mortar.eq.1) then 432 read(15)(xstate(i),i=1,nstate_*mi(1)*(ne+nintpoint)) 433 endif 434 endif 435! 436! face-to-face penalty contact variables 437! 438 if(mortar.eq.1) then 439 read(15) (islavsurf(i),i=1,2*ifacecount+2) 440 read(15) (pslavsurf(i),i=1,3*nintpoint) 441 read(15) (clearini(i),i=1,3*9*ifacecount) 442 endif 443! 444! control parameters 445! 446 read(15) (ctrl(i),i=1,52) 447 read(15) (qaold(i),i=1,2) 448 read(15) output 449 read(15) ttime 450! 451! restart parameters 452! 453 read(15) (irstrt(i),i=1,2) 454! 455 close(15) 456! 457 return 458! 459 15 write(*,*) '*ERROR reading *RESTART,READ: could not open file ', 460 & fnrstrt 461 call exit(201) 462 end 463 464 465 466 467 468 469 470 471 472 473 474 475 476 477 478 479 480 481 482