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 allocation(nload_,nforc_,nboun_,nk_,ne_,nmpc_, 20 & nset_,nalset_,nmat_,ntmat_,npmat_,norien_,nam_,nprint_, 21 & mi,ntrans_,set,meminset,rmeminset,ncs_, 22 & namtot_,ncmat_,memmpc_,ne1d,ne2d,nflow,jobnamec,irstrt, 23 & ithermal,nener,nstate_,irestartstep,inpc,ipoinp,inp, 24 & ntie_,nbody_,nprop_,ipoinpc,nevdamp_,npt_,nslavs,nkon_,mcs, 25 & mortar,ifacecount,nintpoint,infree,nheading_,nobject_, 26 & iuel,iprestr,nstam,ndamp,nef,nbounold,nforcold,nloadold, 27 & nbodyold,mpcend,irobustdesign) 28! 29! calculates a conservative estimate of the size of the 30! fields to be allocated 31! 32! meminset=total # of terms in sets 33! rmeminset=total # of reduced terms (due to use of generate) in 34! sets 35! 36! nstate_ needs only be assigned for 37! a. restart (read from file) 38! b. initial conditions (defined by *depvar) 39! 40 implicit none 41! 42 logical igen,lin,frequency,cyclicsymmetry,composite, 43 & tabular,massflow,beamgeneralsection 44! 45 character*1 selabel,sulabel,inpc(*) 46 character*5 llab 47 character*8 label 48 character*20 mpclabel 49 character*81 set(*),noset,elset,slavset,mastset,noelset,submset, 50 & surface,slavsets,slavsett,mastsets,mastsett,surfset,globalset 51 character*132 jobnamec(*),textpart(16) 52! 53 integer nload_,nforc_,nboun_,nk_,ne_,nmpc_,nset_,nalset_, 54 & nmat_,ntmat_,npmat_,norien_,nam_,nprint_,kode,iline, 55 & istat,n,key,meminset(*),i,js,inoset,mi(*),ii,ipol,inl, 56 & ibounstart,ibounend,ibound,ntrans_,ntmatl,npmatl,ityp,l, 57 & ielset,nope,nteller,nterm,ialset(16),ncs_,rmeminset(*), 58 & islavset,imastset,namtot_,ncmat_,nconstants,memmpc_,j,ipos, 59 & maxrmeminset,ne1d,ne2d,necper,necpsr,necaxr,nesr, 60 & neb32,nn,nflow,nradiate,irestartread,irestartstep,icntrl, 61 & irstrt(*),ithermal(*),nener,nstate_,ipoinp(2,*),inp(3,*), 62 & ntie_,nbody_,nprop_,ipoinpc(0:*),nevdamp_,npt_,nentries, 63 & iposs,iposm,nslavs,nlayer,nkon_,nopeexp,k,iremove,mcs, 64 & ifacecount,nintpoint,mortar,infree(4),nheading_,icfd, 65 & multslav,multmast,nobject_,numnodes,iorientation,id, 66 & irotation,itranslation,nuel,iuel(4,*),number,four, 67 & iprestr,nstam,ier,ndamp,nef,nbounold,nforcold,nloadold, 68 & nbodyold,mpcend,irobustdesign(3),iflag,network,memglobal, 69 & nsubmodel 70! 71 real*8 temperature,tempact,xfreq,tpinc,tpmin,tpmax 72! 73 parameter(nentries=18) 74! 75! icfd=-1: initial value 76! =0: pure mechanical analysis 77! =1: pure CFD analysis 78! =2: mixed mechanical/cfd analysis 79! 80! mi(1): # of integration points 81! mi(2): # of dofs per node 82! mi(3): # of layers in the elements 83! 84 icfd=-1 85! 86 ier=0 87! 88! in the presence of mechanical steps the highest number 89! of DOF is at least 3 90! 91 if(ithermal(2).ne.2) mi(2)=3 92! 93! initialisation of ipoinp 94! 95 do i=1,nentries 96 if(ipoinp(1,i).ne.0) then 97 ipol=i 98 inl=ipoinp(1,i) 99 iline=inp(1,inl)-1 100 exit 101 endif 102 enddo 103! 104 istat=0 105! 106 nsubmodel=0 107 nset_=0 108 maxrmeminset=0 109 necper=0 110 necpsr=0 111 necaxr=0 112 nesr=0 113 neb32=0 114 nradiate=0 115 nkon_=0 116 nuel=0 117! 118 four=4 119! 120 call getnewline(inpc,textpart,istat,n,key,iline,ipol,inl, 121 & ipoinp,inp,ipoinpc) 122 loop: do 123 if(istat.lt.0) then 124 exit 125 endif 126! 127 if(textpart(1)(1:10).eq.'*AMPLITUDE') then 128 nam_=nam_+1 129 do 130 call getnewline(inpc,textpart,istat,n,key,iline,ipol,inl, 131 & ipoinp,inp,ipoinpc) 132 if((istat.lt.0).or.(key.eq.1)) exit 133 namtot_=namtot_+4 134 enddo 135 elseif(textpart(1)(1:19).eq.'*BEAMGENERALSECTION') then 136 mi(3)=max(mi(3),2) 137 do 138 call getnewline(inpc,textpart,istat,n,key,iline,ipol,inl, 139 & ipoinp,inp,ipoinpc) 140 if((istat.lt.0).or.(key.eq.1)) then 141 exit 142 endif 143 nprop_=nprop_+8 144 enddo 145 elseif(textpart(1)(1:12).eq.'*BEAMSECTION') then 146 mi(3)=max(mi(3),2) 147 beamgeneralsection=.false. 148 do i=2,n 149 if((textpart(i)(1:11).eq.'SECTION=BOX').or. 150 & (textpart(i)(1:11).eq.'SECTION=PIP').or. 151 & (textpart(i)(1:11).eq.'SECTION=GEN')) then 152 beamgeneralsection=.true. 153 exit 154 endif 155 enddo 156 if(beamgeneralsection) then 157 do 158 call getnewline(inpc,textpart,istat,n,key,iline,ipol, 159 & inl,ipoinp,inp,ipoinpc) 160 if((istat.lt.0).or.(key.eq.1)) then 161 exit 162 endif 163 nprop_=nprop_+8 164 enddo 165 else 166 call getnewline(inpc,textpart,istat,n,key,iline,ipol, 167 & inl,ipoinp,inp,ipoinpc) 168 endif 169 elseif(textpart(1)(1:10).eq.'*BOUNDARYF') then 170 nam_=nam_+1 171 namtot_=namtot_+1 172 do 173 call getnewline(inpc,textpart,istat,n,key,iline,ipol,inl, 174 & ipoinp,inp,ipoinpc) 175 if((istat.lt.0).or.(key.eq.1)) exit 176! 177 read(textpart(3)(1:10),'(i10)',iostat=istat) ibounstart 178 if(istat.gt.0) then 179 call inputerror(inpc,ipoinpc,iline, 180 & "*BOUNDARYF%",ier) 181 exit 182 endif 183! 184 if(textpart(4)(1:1).eq.' ') then 185 ibounend=ibounstart 186 else 187 read(textpart(4)(1:10),'(i10)',iostat=istat) ibounend 188 if(istat.gt.0) then 189 call inputerror(inpc,ipoinpc,iline, 190 & "*BOUNDARYF%",ier) 191 exit 192 endif 193 endif 194 ibound=ibounend-ibounstart+1 195 ibound=max(1,ibound) 196! 197 read(textpart(1)(1:10),'(i10)',iostat=istat) l 198 if(istat.eq.0) then 199 nboun_=nboun_+ibound 200 if(ntrans_.gt.0) then 201 nmpc_=nmpc_+ibound 202 memmpc_=memmpc_+4*ibound 203 nk_=nk_+1 204 endif 205 else 206 read(textpart(1)(1:80),'(a80)',iostat=istat) elset 207 elset(81:81)=' ' 208 ipos=index(elset,' ') 209! 210! check for element set 211! 212 elset(ipos:ipos)='E' 213c do i=1,nset_ 214c if(set(i).eq.elset) then 215 call cident81(set,elset,nset_,id) 216 i=nset_+1 217 if(id.gt.0) then 218 if(set(id).eq.elset) then 219 i=id 220 nboun_=nboun_+ibound*meminset(i) 221 if(ntrans_.gt.0)then 222 nmpc_=nmpc_+ibound*meminset(i) 223 memmpc_=memmpc_+4*ibound*meminset(i) 224 nk_=nk_+meminset(i) 225 endif 226c exit 227 endif 228 endif 229 if(i.gt.nset_) then 230! 231! check for facial surface 232! 233 elset(ipos:ipos)='T' 234c do i=1,nset_ 235c if(set(i).eq.elset) then 236 call cident81(set,elset,nset_,i) 237 if(i.gt.0) then 238 if(set(i).eq.elset) then 239 nboun_=nboun_+ibound*meminset(i) 240 if(ntrans_.gt.0)then 241 nmpc_=nmpc_+ibound*meminset(i) 242 memmpc_=memmpc_+4*ibound*meminset(i) 243 nk_=nk_+meminset(i) 244 endif 245c exit 246 endif 247 endif 248 endif 249 endif 250 enddo 251 elseif(textpart(1)(1:9).eq.'*BOUNDARY') then 252 nam_=nam_+1 253 namtot_=namtot_+1 254 do 255 call getnewline(inpc,textpart,istat,n,key,iline,ipol,inl, 256 & ipoinp,inp,ipoinpc) 257 if((istat.lt.0).or.(key.eq.1)) exit 258! 259 read(textpart(2)(1:10),'(i10)',iostat=istat) ibounstart 260 if(istat.gt.0) then 261 call inputerror(inpc,ipoinpc,iline, 262 & "*BOUNDARY%",ier) 263 exit 264 endif 265! 266 if(textpart(3)(1:1).eq.' ') then 267 ibounend=ibounstart 268 else 269 read(textpart(3)(1:10),'(i10)',iostat=istat) ibounend 270 if(istat.gt.0) then 271 call inputerror(inpc,ipoinpc,iline, 272 & "*BOUNDARY%",ier) 273 exit 274 endif 275 endif 276 ibound=ibounend-ibounstart+1 277 ibound=max(1,ibound) 278! 279 read(textpart(1)(1:10),'(i10)',iostat=istat) l 280 if(istat.eq.0) then 281 nboun_=nboun_+ibound 282 if(ntrans_.gt.0) then 283 nmpc_=nmpc_+ibound 284 memmpc_=memmpc_+4*ibound 285 nk_=nk_+1 286 endif 287 else 288 read(textpart(1)(1:80),'(a80)',iostat=istat) noset 289 noset(81:81)=' ' 290 ipos=index(noset,' ') 291 noset(ipos:ipos)='N' 292c do i=1,nset_ 293c if(set(i).eq.noset) then 294 call cident81(set,noset,nset_,i) 295 if(i.gt.0) then 296 if(set(i).eq.noset) then 297 nboun_=nboun_+ibound*meminset(i) 298 if(ntrans_.gt.0)then 299 nmpc_=nmpc_+ibound*meminset(i) 300 memmpc_=memmpc_+4*ibound*meminset(i) 301 nk_=nk_+meminset(i) 302 endif 303c exit 304 endif 305 endif 306 endif 307 enddo 308 elseif(textpart(1)(1:4).eq.'*CFD') then 309 iflag=1 310 do i=2,n 311 if(textpart(i)(1:10).eq.'TURBULENCE') then 312 iflag=iflag+1 313 endif 314 enddo 315 if(iflag.eq.2) mi(2)=max(mi(2),6) 316 write(*,*) 'mi(2)= ',mi(2) 317 call getnewline(inpc,textpart,istat,n,key,iline,ipol,inl, 318 & ipoinp,inp,ipoinpc) 319 elseif(textpart(1)(1:6).eq.'*CFLUX') then 320 nam_=nam_+1 321 namtot_=namtot_+1 322 do 323 call getnewline(inpc,textpart,istat,n,key,iline,ipol,inl, 324 & ipoinp,inp,ipoinpc) 325 if((istat.lt.0).or.(key.eq.1)) exit 326! 327 read(textpart(1)(1:10),'(i10)',iostat=istat) l 328 if(istat.eq.0) then 329 nforc_=nforc_+1 330 else 331 read(textpart(1)(1:80),'(a80)',iostat=istat) noset 332 noset(81:81)=' ' 333 ipos=index(noset,' ') 334 noset(ipos:ipos)='N' 335c do i=1,nset_ 336c if(set(i).eq.noset) then 337 call cident81(set,noset,nset_,i) 338 if(i.gt.0) then 339 if(set(i).eq.noset) then 340 nforc_=nforc_+meminset(i) 341c exit 342 endif 343 endif 344 endif 345 enddo 346 elseif(textpart(1)(1:6).eq.'*CLOAD') then 347 nam_=nam_+1 348 namtot_=namtot_+1 349 do 350 call getnewline(inpc,textpart,istat,n,key,iline,ipol,inl, 351 & ipoinp,inp,ipoinpc) 352 if((istat.lt.0).or.(key.eq.1)) exit 353! 354 read(textpart(1)(1:10),'(i10)',iostat=istat) l 355 if(istat.eq.0) then 356 if(ntrans_.eq.0) then 357 nforc_=nforc_+1 358 else 359 nforc_=nforc_+3 360 endif 361 else 362 read(textpart(1)(1:80),'(a80)',iostat=istat) noset 363 noset(81:81)=' ' 364 ipos=index(noset,' ') 365 noset(ipos:ipos)='N' 366c do i=1,nset_ 367c if(set(i).eq.noset) then 368 call cident81(set,noset,nset_,i) 369 if(i.gt.0) then 370 if(set(i).eq.noset) then 371 if(ntrans_.eq.0) then 372 nforc_=nforc_+meminset(i) 373 else 374 nforc_=nforc_+3*meminset(i) 375 endif 376c exit 377 endif 378 endif 379 endif 380 enddo 381 elseif((textpart(1)(1:13).eq.'*CONDUCTIVITY').or. 382 & (textpart(1)(1:8).eq.'*DENSITY').or. 383 & (textpart(1)(1:10).eq.'*EXPANSION').or. 384 & (textpart(1)(1:15).eq.'*FLUIDCONSTANTS').or. 385 & (textpart(1)(1:13).eq.'*SPECIFICHEAT').or. 386 & (textpart(1)(1:23).eq.'*ELECTRICALCONDUCTIVITY')) then 387 ntmatl=0 388 do 389 call getnewline(inpc,textpart,istat,n,key,iline,ipol,inl, 390 & ipoinp,inp,ipoinpc) 391 if((istat.lt.0).or.(key.eq.1)) exit 392 ntmatl=ntmatl+1 393 ntmat_=max(ntmatl,ntmat_) 394 enddo 395 elseif(textpart(1)(1:11).eq.'*CONSTRAINT') then 396 do 397 call getnewline(inpc,textpart,istat,n,key,iline,ipol,inl, 398 & ipoinp,inp,ipoinpc) 399 if((istat.lt.0).or.(key.eq.1)) exit 400 nobject_=nobject_+1 401 enddo 402 elseif(textpart(1)(1:15).eq.'*CONTACTDAMPING') then 403 ncmat_=max(8,ncmat_) 404 ntmat_=max(1,ntmat_) 405 call getnewline(inpc,textpart,istat,n,key,iline,ipol,inl, 406 & ipoinp,inp,ipoinpc) 407 elseif(textpart(1)(1:12).eq.'*CONTACTPAIR') then 408 do 409 call getnewline(inpc,textpart,istat,n,key,iline,ipol,inl, 410 & ipoinp,inp,ipoinpc) 411 if((istat.lt.0).or.(key.eq.1)) exit 412 ntie_=ntie_+1 413 enddo 414 elseif(textpart(1)(1:13).eq.'*CONTACTPRINT') then 415 do 416 call getnewline(inpc,textpart,istat,n,key,iline,ipol,inl, 417 & ipoinp,inp,ipoinpc) 418 if((istat.lt.0).or.(key.eq.1)) exit 419 nprint_=nprint_+n 420 enddo 421 elseif(textpart(1)(1:9).eq.'*COUPLING') then 422 surface(1:1)=' ' 423 iorientation=0 424 do i=2,n 425 if(textpart(i)(1:8).eq.'SURFACE=') then 426 surface=textpart(i)(9:88) 427 ipos=index(surface,' ') 428 surface(ipos:ipos)='T' 429 elseif(textpart(i)(1:12).eq.'ORIENTATION=') then 430 iorientation=1 431 endif 432 enddo 433 if(surface(1:1).ne.' ') then 434c do i=1,nset_ 435c surface(ipos:ipos)='T' 436c if(set(i).eq.surface) then 437c numnodes=8*meminset(i) 438c exit 439c endif 440c surface(ipos:ipos)='S' 441c if(set(i).eq.surface) then 442c numnodes=meminset(i) 443c exit 444c endif 445c enddo 446 surface(ipos:ipos)='T' 447 call cident81(set,surface,nset_,i) 448 if(i.gt.0) then 449 if(set(i).eq.surface) then 450 numnodes=8*meminset(i) 451 endif 452 endif 453 surface(ipos:ipos)='S' 454 call cident81(set,surface,nset_,i) 455 if(i.gt.0) then 456 if(set(i).eq.surface) then 457 numnodes=meminset(i) 458 endif 459 endif 460 endif 461 call getnewline(inpc,textpart,istat,n,key,iline,ipol,inl, 462 & ipoinp,inp,ipoinpc) 463 elseif(textpart(1)(1:6).eq.'*CREEP') then 464 ntmatl=0 465 npmat_=max(2,npmat_) 466 if(ncmat_.le.2) then 467! elastic isotropic 468 ncmat_=max(9,ncmat_) 469 else 470! elastic anisotropic 471 ncmat_=max(19,ncmat_) 472 endif 473 do 474 call getnewline(inpc,textpart,istat,n,key,iline,ipol,inl, 475 & ipoinp,inp,ipoinpc) 476 if((istat.lt.0).or.(key.eq.1)) exit 477 ntmatl=ntmatl+1 478 enddo 479 ntmat_=max(ntmatl,ntmat_) 480 elseif(textpart(1)(1:16).eq.'*CYCLICHARDENING') then 481 ntmatl=0 482 do 483 call getnewline(inpc,textpart,istat,n,key,iline,ipol,inl, 484 & ipoinp,inp,ipoinpc) 485 if((istat.lt.0).or.(key.eq.1)) exit 486 read(textpart(3)(1:20),'(f20.0)',iostat=istat) 487 & temperature 488 if(istat.gt.0) then 489 call inputerror(inpc,ipoinpc,iline, 490 & "*CYCLIC HARDENING%",ier) 491 exit 492 endif 493 if(ntmatl.eq.0) then 494 npmatl=0 495 ntmatl=ntmatl+1 496 ntmat_=max(ntmatl,ntmat_) 497 tempact=temperature 498 elseif(temperature.ne.tempact) then 499 npmatl=0 500 ntmatl=ntmatl+1 501 ntmat_=max(ntmatl,ntmat_) 502 tempact=temperature 503 endif 504 npmatl=npmatl+1 505 npmat_=max(npmatl,npmat_) 506 enddo 507 elseif(textpart(1)(1:20).eq.'*CYCLICSYMMETRYMODEL') then 508! 509! possible MPC's: static temperature, displacements(velocities) 510! and static pressure 511! 512 nk_=nk_+1 513c nmpc_=nmpc_+5*ncs_ 514c memmpc_=memmpc_+125*ncs_ 515! 516! change on 11th of Dec. 2020 517! 518 nmpc_=nmpc_+(mi(2)+1)*ncs_ 519 memmpc_=memmpc_+25*(mi(2)+1)*ncs_ 520 ntrans_=ntrans_+1 521 do 522 call getnewline(inpc,textpart,istat,n,key,iline,ipol,inl, 523 & ipoinp,inp,ipoinpc) 524 if((istat.lt.0).or.(key.eq.1)) exit 525 enddo 526 elseif(textpart(1)(1:8).eq.'*DAMPING') then 527 do i=2,n 528 if(textpart(i)(1:11).eq.'STRUCTURAL=') then 529 ndamp=1 530 exit 531 endif 532 enddo 533 call getnewline(inpc,textpart,istat,n,key,iline,ipol,inl, 534 & ipoinp,inp,ipoinpc) 535 elseif(textpart(1)(1:8).eq.'*DASHPOT') then 536 nmat_=nmat_+1 537 frequency=.false. 538 call getnewline(inpc,textpart,istat,n,key,iline,ipol,inl, 539 & ipoinp,inp,ipoinpc) 540 if((istat.lt.0).or.(key.eq.1)) then 541 call inputerror(inpc,ipoinpc,iline, 542 & "*DASHPOT%",ier) 543 cycle 544 endif 545 read(textpart(2)(1:20),'(f20.0)',iostat=istat) 546 & xfreq 547 if(istat.gt.0) then 548 call inputerror(inpc,ipoinpc,iline, 549 & "*DASHPOT%",ier) 550 cycle 551 endif 552 if(xfreq.gt.0.d0) frequency=.true. 553 iline=iline-1 554 if(.not.frequency) then 555 ntmatl=0 556 ncmat_=max(2,ncmat_) 557 do 558 call getnewline(inpc,textpart,istat,n,key,iline,ipol, 559 & inl,ipoinp,inp,ipoinpc) 560 if((istat.lt.0).or.(key.eq.1)) exit 561 ntmatl=ntmatl+1 562 ntmat_=max(ntmatl,ntmat_) 563 enddo 564 else 565 ntmatl=0 566 do 567 call getnewline(inpc,textpart,istat,n,key,iline,ipol, 568 & inl,ipoinp,inp,ipoinpc) 569 if((istat.lt.0).or.(key.eq.1)) exit 570 read(textpart(3)(1:20),'(f20.0)',iostat=istat) 571 & temperature 572 if(istat.gt.0) then 573 call inputerror(inpc,ipoinpc,iline, 574 & "*DASHPOT%",ier) 575 exit 576 endif 577 if(ntmatl.eq.0) then 578 npmatl=0 579 ntmatl=ntmatl+1 580 ntmat_=max(ntmatl,ntmat_) 581 tempact=temperature 582 elseif(temperature.ne.tempact) then 583 npmatl=0 584 ntmatl=ntmatl+1 585 ntmat_=max(ntmatl,ntmat_) 586 tempact=temperature 587 endif 588 npmatl=npmatl+1 589 npmat_=max(npmatl,npmat_) 590 enddo 591 if(ncmat_.ge.9) ncmat_=max(19,ncmat_) 592 endif 593 elseif(textpart(1)(1:22).eq.'*DEFORMATIONPLASTICITY') then 594 ncmat_=max(5,ncmat_) 595 ntmatl=0 596 do 597 call getnewline(inpc,textpart,istat,n,key,iline,ipol,inl, 598 & ipoinp,inp,ipoinpc) 599 if((istat.lt.0).or.(key.eq.1)) exit 600 ntmatl=ntmatl+1 601 ntmat_=max(ntmatl,ntmat_) 602 enddo 603 elseif(textpart(1)(1:7).eq.'*DEPVAR') then 604 do 605 call getnewline(inpc,textpart,istat,n,key,iline,ipol,inl, 606 & ipoinp,inp,ipoinpc) 607 if((istat.lt.0).or.(key.eq.1)) exit 608 read(textpart(1)(1:10),'(i10)',iostat=istat) l 609 if(istat.lt.0) exit 610 nstate_=max(l,nstate_) 611 enddo 612 elseif(textpart(1)(1:16).eq.'*DESIGNVARIABLES') then 613 ntie_=ntie_+1 614 call getnewline(inpc,textpart,istat,n,key,iline,ipol,inl, 615 & ipoinp,inp,ipoinpc) 616 elseif(textpart(1)(1:15).eq.'*DESIGNRESPONSE') then 617 do 618 call getnewline(inpc,textpart,istat,n,key,iline,ipol,inl, 619 & ipoinp,inp,ipoinpc) 620 if((istat.lt.0).or.(key.eq.1)) exit 621 nobject_=nobject_+1 622 enddo 623 elseif(textpart(1)(1:21).eq.'*DISTRIBUTINGCOUPLING') then 624 nmpc_=nmpc_+3 625 memmpc_=memmpc_+3 626 do 627 call getnewline(inpc,textpart,istat,n,key,iline,ipol,inl, 628 & ipoinp,inp,ipoinpc) 629 if((istat.lt.0).or.(key.eq.1)) exit 630! 631 read(textpart(1)(1:10),'(i10)',iostat=istat) l 632 if(istat.eq.0) then 633 memmpc_=memmpc_+3 634 else 635 read(textpart(1)(1:80),'(a80)',iostat=istat) noset 636 noset(81:81)=' ' 637 ipos=index(noset,' ') 638 noset(ipos:ipos)='N' 639c do i=1,nset_ 640c if(set(i).eq.noset) then 641 call cident81(set,noset,nset_,i) 642 if(i.gt.0) then 643 if(set(i).eq.noset) then 644 memmpc_=memmpc_+3*meminset(i) 645c exit 646 endif 647 endif 648 endif 649 enddo 650 elseif(textpart(1)(2:13).eq.'DISTRIBUTING') then 651 irotation=0 652 itranslation=0 653 do 654 call getnewline(inpc,textpart,istat,n,key,iline,ipol,inl, 655 & ipoinp,inp,ipoinpc) 656 if((istat.lt.0).or.(key.eq.1)) exit 657! 658 read(textpart(1)(1:10),'(i10)',iostat=istat) ibounstart 659 if(istat.gt.0) then 660 call inputerror(inpc,ipoinpc,iline, 661 & "*BOUNDARY%",ier) 662 exit 663 endif 664! 665 if(textpart(2)(1:1).eq.' ') then 666 ibounend=ibounstart 667 else 668 read(textpart(2)(1:10),'(i10)',iostat=istat) ibounend 669 if(istat.gt.0) then 670 call inputerror(inpc,ipoinpc,iline, 671 & "*BOUNDARY%",ier) 672 exit 673 endif 674 endif 675 ibounstart=max(4,ibounstart) 676 ibounend=min(6,ibounend) 677 ibound=max(0,ibounend-ibounstart+1) 678! 679 if(itranslation.eq.0) then 680! 681! translational dofs 3 MPC's + a two-term MPC for each 682! participating node 683! 684 npt_=max(npt_,numnodes) 685! 686 nmpc_=nmpc_+3*npt_+3 687 memmpc_=memmpc_+6*npt_+3*(npt_+1) 688 nk_=nk_+npt_ 689 itranslation=1 690 endif 691! 692! rotational dofs 693! 694 if(ibound.gt.0) then 695 if(irotation.eq.0) then 696! 697! a MPC connecting the dofs 4-6 to dofs 1-3 of 698! a rotational node; generation of a inhomogeneous 699! node 700! 701 nmpc_=nmpc_+3 702 memmpc_=memmpc_+6 703 nk_=nk_+4 704 irotation=1 705 endif 706 nmpc_=nmpc_+ibound 707 memmpc_=memmpc_+ibound*(3*npt_+2) 708 nboun_=nboun_+ibound 709 endif 710 enddo 711! 712 elseif(textpart(1)(2:18).eq.'DISTRIBUTIONTABLE') then 713 call getnewline(inpc,textpart,istat,n,key,iline,ipol,inl, 714 & ipoinp,inp,ipoinpc) 715 elseif(textpart(1)(2:13).eq.'DISTRIBUTION') then 716 do 717 call getnewline(inpc,textpart,istat,n,key,iline,ipol,inl, 718 & ipoinp,inp,ipoinpc) 719 if((istat.lt.0).or.(key.eq.1)) exit 720 norien_=norien_+1 721 enddo 722 elseif((textpart(1)(1:6).eq.'*DLOAD').or. 723 & (textpart(1)(1:7).eq.'*DSLOAD').or. 724 & (textpart(1)(1:6).eq.'*DFLUX').or. 725 & (textpart(1)(1:9).eq.'*MASSFLOW').or. 726 & (textpart(1)(1:5).eq.'*FILM')) then 727 massflow=.false. 728 if((textpart(1)(1:5).ne.'*FILM').and. 729 & (textpart(1)(1:9).ne.'*MASSFLOW')) then 730 nam_=nam_+1 731 namtot_=namtot_+1 732 elseif(textpart(1)(1:9).ne.'*MASSFLOW') then 733 nam_=nam_+2 734 namtot_=namtot_+2 735 else 736 massflow=.true. 737 endif 738 do 739 call getnewline(inpc,textpart,istat,n,key,iline,ipol,inl, 740 & ipoinp,inp,ipoinpc) 741 if((istat.lt.0).or.(key.eq.1)) exit 742 read(textpart(2)(1:5),'(a5)',iostat=istat) llab 743 if((llab.eq.'GRAV ').or.(llab.eq.'CENTR').or. 744 & (llab.eq.'NEWTO')) then 745 nbody_=nbody_+1 746 cycle 747 endif 748 read(textpart(1)(1:10),'(i10)',iostat=istat) l 749 if(istat.eq.0) then 750 nload_=nload_+1 751 if(massflow) then 752 nmpc_=nmpc_+1 753 memmpc_=memmpc_+3 754 endif 755 else 756 read(textpart(1)(1:80),'(a80)',iostat=istat) elset 757 elset(81:81)=' ' 758 ipos=index(elset,' ') 759! 760! check for element set 761! 762 elset(ipos:ipos)='E' 763c do i=1,nset_ 764c if(set(i).eq.elset) then 765 call cident81(set,elset,nset_,id) 766 i=nset_+1 767 if(id.gt.0) then 768 if(set(id).eq.elset) then 769 i=id 770 nload_=nload_+meminset(i) 771 if(massflow) then 772 nmpc_=nmpc_+meminset(i) 773 memmpc_=memmpc_+3*meminset(i) 774 endif 775c exit 776 endif 777 endif 778 if(i.gt.nset_) then 779! 780! check for facial surface 781! 782 elset(ipos:ipos)='T' 783c do i=1,nset_ 784c if(set(i).eq.elset) then 785 call cident81(set,elset,nset_,i) 786 if(i.gt.0) then 787 if(set(i).eq.elset) then 788 nload_=nload_+meminset(i) 789 if(massflow) then 790 nmpc_=nmpc_+meminset(i) 791 memmpc_=memmpc_+3*meminset(i) 792 endif 793c exit 794 endif 795 endif 796 endif 797 endif 798 enddo 799 elseif((textpart(1)(1:8).eq.'*DYNAMIC').or. 800 & (textpart(1)(1:32).eq.'*COUPLEDTEMPERATURE-DISPLACEMENT') 801 & .or. 802 & (textpart(1)(1:34).eq. 803 & '*UNCOUPLEDTEMPERATURE-DISPLACEMENT'))then 804! 805! change of number of integration points except for a pure 806! CFD-calculation 807! 808 if(icfd.ne.1) then 809 if((mi(1).eq.1).or.(mi(1).eq.8)) then 810 mi(1)=27 811 elseif(mi(1).eq.4) then 812 mi(1)=15 813 elseif(mi(1).eq.2) then 814 mi(1)=9 815 endif 816 endif 817 call getnewline(inpc,textpart,istat,n,key,iline,ipol,inl, 818 & ipoinp,inp,ipoinpc) 819 elseif(textpart(1)(1:8).eq.'*ELPRINT') then 820 do 821 call getnewline(inpc,textpart,istat,n,key,iline,ipol,inl, 822 & ipoinp,inp,ipoinpc) 823 if((istat.lt.0).or.(key.eq.1)) exit 824 nprint_=nprint_+n 825 enddo 826 elseif(textpart(1)(1:8).eq.'*ELASTIC') then 827 ntmatl=0 828 ityp=2 829 ncmat_=max(2,ncmat_) 830 do i=2,n 831 if(textpart(i)(1:5).eq.'TYPE=') then 832 if(textpart(i)(6:8).eq.'ISO') then 833 ityp=2 834 ncmat_=max(2,ncmat_) 835 elseif((textpart(i)(6:10).eq.'ORTHO').or. 836 & (textpart(i)(6:10).eq.'ENGIN')) then 837 ityp=9 838 ncmat_=max(9,ncmat_) 839 elseif(textpart(i)(6:10).eq.'ANISO') then 840 ityp=21 841 ncmat_=max(21,ncmat_) 842 endif 843 exit 844 endif 845 enddo 846 if(ityp.eq.2) then 847 do 848 call getnewline(inpc,textpart,istat,n,key,iline,ipol, 849 & inl,ipoinp,inp,ipoinpc) 850 if((istat.lt.0).or.(key.eq.1)) exit 851 ntmatl=ntmatl+1 852 enddo 853 ntmat_=max(ntmatl,ntmat_) 854 elseif(ityp.eq.9) then 855 do 856 call getnewline(inpc,textpart,istat,n,key,iline,ipol, 857 & inl,ipoinp,inp,ipoinpc) 858 if((istat.lt.0).or.(key.eq.1)) exit 859 ntmatl=ntmatl+1 860 iline=iline+1 861 enddo 862 ntmat_=max(ntmatl,ntmat_) 863 elseif(ityp.eq.21) then 864 do 865 call getnewline(inpc,textpart,istat,n,key,iline,ipol, 866 & inl,ipoinp,inp,ipoinpc) 867 if((istat.lt.0).or.(key.eq.1)) exit 868 ntmatl=ntmatl+1 869 iline=iline+2 870 enddo 871 ntmat_=max(ntmatl,ntmat_) 872 endif 873 elseif(textpart(1)(1:17).eq.'*ELECTROMAGNETICS') then 874 mi(2)=max(mi(2),5) 875 call getnewline(inpc,textpart,istat,n,key,iline,ipol,inl, 876 & ipoinp,inp,ipoinpc) 877 elseif((textpart(1)(1:8).eq.'*ELEMENT').and. 878 & (textpart(1)(1:14).ne.'*ELEMENTOUTPUT')) then 879 ielset=0 880! 881 loop1: do i=2,n 882 if(textpart(i)(1:6).eq.'ELSET=') then 883 elset=textpart(i)(7:86) 884 elset(81:81)=' ' 885 ipos=index(elset,' ') 886 elset(ipos:ipos)='E' 887 ielset=1 888c do js=1,nset_ 889c if(set(js).eq.elset) exit 890c enddo 891c if(js.gt.nset_) then 892c nset_=nset_+1 893c set(nset_)=elset 894c endif 895 call cident81(set,elset,nset_,id) 896 js=nset_+1 897 if(id.gt.0) then 898 if(set(id).eq.elset) js=id 899 endif 900 if(js.gt.nset_) then 901 nset_=nset_+1 902 do j=nset_,id+2,-1 903 meminset(j)=meminset(j-1) 904 rmeminset(j)=rmeminset(j-1) 905 set(j)=set(j-1) 906 enddo 907 js=id+1 908 set(js)=elset 909 meminset(js)=0 910 rmeminset(js)=0 911 endif 912 elseif(textpart(i)(1:5).eq.'TYPE=') then 913 read(textpart(i)(6:13),'(a8)') label 914 if(label.eq.' ') then 915 write(*,*) 916 & '*ERROR in allocation: element type is lacking' 917 write(*,*) ' ' 918 call inputerror(inpc,ipoinpc,iline, 919 & "*ELEMENT or *ELEMENT OUTPUT%",ier) 920 exit 921 endif 922 if((label(1:2).eq.'DC').and.(label(1:7).ne.'DCOUP3D')) 923 & then 924 label(1:7)=label(2:8) 925 label(8:8)=' ' 926 endif 927! 928 nopeexp=0 929! 930 if(label.eq.'C3D20 ') then 931 mi(1)=max(mi(1),27) 932 nope=20 933 nopeexp=20 934 elseif(label(1:8).eq.'C3D20R ') then 935 mi(1)=max(mi(1),8) 936 nope=20 937 nopeexp=20 938 elseif((label.eq.'C3D8R ').or.(label.eq.'F3D8R ')) 939 & then 940 mi(1)=max(mi(1),1) 941 nope=8 942 nopeexp=8 943 elseif((label.eq.'C3D10 ').or. 944 & (label.eq.'C3D10T ')) then 945 mi(1)=max(mi(1),4) 946 nope=10 947 nopeexp=10 948 elseif((label.eq.'C3D4 ').or. 949 & (label.eq.'F3D4 ')) then 950 mi(1)=max(mi(1),1) 951 nope=4 952 nopeexp=4 953 elseif(label.eq.'C3D15 ') then 954 mi(1)=max(mi(1),9) 955 nope=15 956 nopeexp=15 957 elseif(label.eq.'C3D6 ') then 958 mi(1)=max(mi(1),2) 959 nope=6 960 nopeexp=6 961 elseif(label.eq.'F3D6 ') then 962 mi(1)=max(mi(1),1) 963 nope=6 964 nopeexp=6 965 elseif((label.eq.'C3D8 ').or.(label.eq.'F3D8 ')) 966 & then 967 mi(1)=max(mi(1),8) 968 nope=8 969 nopeexp=8 970c Bernhardi start 971 elseif(label.eq.'C3D8I ') then 972 mi(1)=max(mi(1),8) 973 nope=8 974 nopeexp=11 975c Bernhardi end 976 elseif((label.eq.'CPE3 ').or. 977 & (label.eq.'CPS3 ').or. 978 & (label.eq.'CAX3 ').or. 979 & (label.eq.'M3D3 ').or. 980 & (label.eq.'S3 ')) then 981 mi(1)=max(mi(1),2) 982 nope=3 983 nopeexp=9 984 elseif((label.eq.'CPE4R ').or. 985 & (label.eq.'CPS4R ').or. 986 & (label.eq.'CAX4R ').or. 987 & (label.eq.'M3D4R ').or. 988 & (label.eq.'S4R ')) then 989 mi(1)=max(mi(1),1) 990 nope=4 991 nopeexp=12 992 elseif((label.eq.'CPE4 ').or. 993 & (label.eq.'CPS4 ').or. 994 & (label.eq.'CAX4 ').or. 995 & (label.eq.'M3D4 ')) then 996 mi(1)=max(mi(1),8) 997 nope=4 998 nopeexp=12 999 elseif(label.eq.'S4 ') then 1000 mi(1)=max(mi(1),8) 1001 nope=4 1002! modified into C3D8I (11 nodes) 1003 nopeexp=15 1004 elseif((label.eq.'CPE6 ').or. 1005 & (label.eq.'CPS6 ').or. 1006 & (label.eq.'CAX6 ').or. 1007 & (label.eq.'M3D6 ').or. 1008 & (label.eq.'S6 ')) then 1009 mi(1)=max(mi(1),9) 1010 nope=6 1011 nopeexp=21 1012 elseif((label.eq.'CPE8R ').or. 1013 & (label.eq.'CPS8R ').or. 1014 & (label.eq.'CAX8R ').or. 1015 & (label.eq.'M3D8R ').or. 1016 & (label.eq.'S8R ')) then 1017 mi(1)=max(mi(1),8) 1018 nope=8 1019 nopeexp=28 1020 elseif((label.eq.'CPE8 ').or. 1021 & (label.eq.'CPS8 ').or. 1022 & (label.eq.'CAX8 ').or. 1023 & (label.eq.'M3D8 ').or. 1024 & (label.eq.'S8 ')) then 1025 mi(1)=max(mi(1),27) 1026 nope=8 1027 nopeexp=28 1028 elseif((label.eq.'B31 ').or. 1029 & (label.eq.'B21 ').or. 1030 & (label.eq.'T3D2 ').or. 1031 & (label.eq.'T2D2 ')) then 1032 mi(1)=max(mi(1),8) 1033 mi(3)=max(mi(3),2) 1034 nope=2 1035! modified into C3D8I (11 nodes) 1036 nopeexp=13 1037 elseif(label.eq.'B31R ') then 1038 mi(1)=max(mi(1),1) 1039 nope=2 1040 nopeexp=10 1041 elseif((label.eq.'B32 ').or. 1042 & (label.eq.'T3D3 ')) then 1043 mi(1)=max(mi(1),27) 1044 mi(3)=max(mi(3),2) 1045 nope=3 1046 nopeexp=23 1047 elseif(label.eq.'B32R ') then 1048 mi(1)=max(mi(1),50) 1049 nope=3 1050 nopeexp=23 1051 elseif(label(1:8).eq.'DASHPOTA') then 1052 label='EDSHPTA1' 1053 nope=2 1054 nopeexp=2 1055 elseif(label(1:7).eq.'DCOUP3D') then 1056 nope=1 1057 nopeexp=1 1058 elseif(label(1:1).eq.'D') then 1059 nope=3 1060 nopeexp=3 1061 mi(2)=max(3,mi(2)) 1062 elseif(label(1:7).eq.'SPRINGA') then 1063 mi(1)=max(mi(1),1) 1064 label='ESPRNGA1' 1065 nope=2 1066 nopeexp=2 1067 elseif(label(1:7).eq.'SPRING1') then 1068 mi(1)=max(mi(1),1) 1069 label='ESPRNG10' 1070 nope=1 1071 nopeexp=1 1072 ncmat_=max(3,ncmat_) 1073 elseif(label(1:7).eq.'SPRING2') then 1074 mi(1)=max(mi(1),1) 1075 label='ESPRNG21' 1076 nope=2 1077 nopeexp=2 1078 ncmat_=max(4,ncmat_) 1079 elseif(label.eq.'GAPUNI ') then 1080 mi(1)=max(mi(1),1) 1081 label='ESPGAPA1' 1082 nope=2 1083 nopeexp=2 1084 elseif(label(1:4).eq.'MASS') then 1085 nope=1 1086 nopeexp=1 1087 elseif(label(1:1).eq.'U') then 1088! 1089! the number uniquely characterizes the 1090! element name (consisting of 4 freely 1091! selectable characters in position 2..5) 1092! 1093 number=ichar(label(2:2))*256**3+ 1094 & ichar(label(3:3))*256**2+ 1095 & ichar(label(4:4))*256+ 1096 & ichar(label(5:5)) 1097 nope=-1 1098 call nidentk(iuel,number,nuel,id,four) 1099 if(id.gt.0) then 1100 if(iuel(1,id).eq.number) then 1101 mi(1)=max(mi(1),iuel(2,id)) 1102 mi(2)=max(mi(2),iuel(3,id)) 1103 nope=iuel(4,id) 1104 nopeexp=nope 1105 endif 1106 endif 1107 if(nope.eq.-1) then 1108 write(*,*) '*ERROR reading *ELEMENT' 1109 write(*,*) ' nonexistent element type:' 1110 write(*,*) ' ',label 1111 call inputerror(inpc,ipoinpc,iline, 1112 & "*ELEMENT%",ier) 1113 cycle loop 1114 endif 1115 endif 1116 if(label(1:1).eq.'F') then 1117 mi(2)=max(mi(2),4) 1118 if(icfd.eq.-1) then 1119 icfd=1 1120 elseif(icfd.eq.0) then 1121 icfd=2 1122 endif 1123 else 1124 if(icfd.eq.-1) then 1125 icfd=0 1126 elseif(icfd.eq.1) then 1127 icfd=2 1128 endif 1129 endif 1130 endif 1131 enddo loop1 1132! 1133 loop2:do 1134 call getnewline(inpc,textpart,istat,n,key,iline,ipol,inl, 1135 & ipoinp,inp,ipoinpc) 1136 if((istat.lt.0).or.(key.eq.1)) exit 1137 read(textpart(1)(1:10),'(i10)',iostat=istat) i 1138 if(istat.gt.0) then 1139 call inputerror(inpc,ipoinpc,iline, 1140 & "*ELEMENT or *ELEMENT OUTPUT%",ier) 1141 exit 1142 endif 1143c Bernhardi start 1144c space for incompatible mode nodes 1145 if(label(1:5).eq.'C3D8I') then 1146 nk_=nk_+3 1147 endif 1148c Bernhardi end 1149 if(label(1:2).ne.'C3') then 1150 if(label(1:3).eq.'CPE') then 1151 necper=necper+1 1152 elseif(label(1:2).eq.'CP') then 1153 necpsr=necpsr+1 1154 elseif(label(1:1).eq.'C') then 1155 necaxr=necaxr+1 1156 elseif((label(1:1).eq.'S').or. 1157 & ((label(1:1).eq.'M').and.(label(1:4).ne.'MASS'))) 1158 & then 1159 nesr=nesr+1 1160 elseif((label(1:1).eq.'B').or. 1161 & (label(1:1).eq.'T')) then 1162 neb32=neb32+1 1163 elseif(label(1:1).eq.'D') then 1164 nflow=nflow+1 1165 elseif(label(1:1).eq.'F') then 1166 nef=nef+1 1167 endif 1168 endif 1169 nteller=n-1 1170 if(nteller.lt.nope) then 1171 do 1172 call getnewline(inpc,textpart,istat,n,key,iline, 1173 & ipol,inl,ipoinp,inp,ipoinpc) 1174 if((istat.lt.0).or.(key.eq.1)) exit loop2 1175 if(nteller+n.gt.nope) n=nope-nteller 1176 nteller=nteller+n 1177 if(nteller.eq.nope) exit 1178 enddo 1179 endif 1180 ne_=max(ne_,i) 1181 nkon_=nkon_+nopeexp 1182 if(ielset.eq.1) then 1183 meminset(js)=meminset(js)+1 1184 rmeminset(js)=rmeminset(js)+1 1185 endif 1186c ! 1187c ! up to 8 new mpc's with 22 terms in each mpc 1188c ! (21 = 7 nodes x 3 dofs + inhomogeneous term) 1189c ! 1190 enddo loop2 1191 elseif((textpart(1)(1:5).eq.'*NSET').or. 1192 & (textpart(1)(1:6).eq.'*ELSET')) then 1193 if(textpart(1)(1:5).eq.'*NSET') 1194 & then 1195 noelset=textpart(2)(6:85) 1196 noelset(81:81)=' ' 1197 ipos=index(noelset,' ') 1198 noelset(ipos:ipos)='N' 1199 kode=0 1200 else 1201 noelset=textpart(2)(7:86) 1202 noelset(81:81)=' ' 1203 ipos=index(noelset,' ') 1204 noelset(ipos:ipos)='E' 1205 kode=1 1206 endif 1207! 1208! check whether new set name or old one 1209! 1210c do js=1,nset_ 1211c if(set(js).eq.noelset) exit 1212c enddo 1213c if(js.gt.nset_) then 1214c nset_=nset_+1 1215c set(nset_)=noelset 1216c nn=nset_ 1217c else 1218c nn=js 1219c endif 1220 call cident81(set,noelset,nset_,id) 1221 nn=nset_+1 1222 if(id.gt.0) then 1223 if(set(id).eq.noelset) nn=id 1224 endif 1225 if(nn.gt.nset_) then 1226 nset_=nset_+1 1227 do j=nset_,id+2,-1 1228 meminset(j)=meminset(j-1) 1229 rmeminset(j)=rmeminset(j-1) 1230 set(j)=set(j-1) 1231 enddo 1232 nn=id+1 1233 set(nn)=noelset 1234 meminset(nn)=0 1235 rmeminset(nn)=0 1236 endif 1237! 1238 if((n.gt.2).and.(textpart(3)(1:8).eq.'GENERATE')) then 1239 igen=.true. 1240 else 1241 igen=.false. 1242 endif 1243 do 1244 call getnewline(inpc,textpart,istat,n,key,iline,ipol,inl, 1245 & ipoinp,inp,ipoinpc) 1246 if((istat.lt.0).or.(key.eq.1)) exit 1247 if(igen) then 1248 if(textpart(2)(1:1).eq.' ') 1249 & textpart(2)=textpart(1) 1250 if(textpart(3)(1:1).eq.' ') 1251 & textpart(3)='1 ' 1252 do i=1,3 1253 read(textpart(i)(1:10),'(i10)',iostat=istat) 1254 & ialset(i) 1255 if(istat.gt.0) then 1256 call inputerror(inpc,ipoinpc,iline, 1257 & "*NSET or *ELSET%",ier) 1258 exit 1259 endif 1260 enddo 1261 meminset(nn)=meminset(nn)+ 1262 & (ialset(2)-ialset(1))/ialset(3)+1 1263 rmeminset(nn)=rmeminset(nn)+3 1264 else 1265 do i=1,n 1266 read(textpart(i)(1:10),'(i10)',iostat=istat) 1267 & ialset(i) 1268 if(istat.gt.0) then 1269 noelset=textpart(i)(1:80) 1270 noelset(81:81)=' ' 1271 ipos=index(noelset,' ') 1272 if(kode.eq.0) then 1273 noelset(ipos:ipos)='N' 1274 else 1275 noelset(ipos:ipos)='E' 1276 endif 1277c do j=1,nset_ 1278c if(noelset.eq.set(j)) then 1279 call cident81(set,noelset,nset_,j) 1280 if(j.gt.0) then 1281 if(set(j).eq.noelset) then 1282 meminset(nn)=meminset(nn)+ 1283 & meminset(j) 1284 rmeminset(nn)=rmeminset(nn)+ 1285 & rmeminset(j) 1286c exit 1287 endif 1288 endif 1289 else 1290 meminset(nn)=meminset(nn)+1 1291 rmeminset(nn)=rmeminset(nn)+1 1292 endif 1293 enddo 1294 endif 1295 enddo 1296 elseif((textpart(1)(1:9).eq.'*EQUATION').or. 1297 & (textpart(1)(1:10).eq.'*EQUATIONF')) then 1298 iremove=0 1299 do i=2,n 1300 if(textpart(i)(1:6).eq.'REMOVE') iremove=1 1301 enddo 1302 do 1303 call getnewline(inpc,textpart,istat,n,key,iline,ipol,inl, 1304 & ipoinp,inp,ipoinpc) 1305 if(iremove.eq.1) exit 1306 if((istat.lt.0).or.(key.eq.1)) exit 1307 read(textpart(1)(1:10),'(i10)',iostat=istat) nterm 1308 if(ntrans_.eq.0) then 1309 nmpc_=nmpc_+1 1310 memmpc_=memmpc_+nterm 1311 else 1312 nmpc_=nmpc_+3 1313 memmpc_=memmpc_+3*nterm 1314 endif 1315 ii=0 1316 do 1317 call getnewline(inpc,textpart,istat,n,key,iline,ipol, 1318 & inl,ipoinp,inp,ipoinpc) 1319 if((istat.lt.0).or.(key.eq.1)) exit 1320 ii=ii+n/3 1321 if(ii.eq.nterm) exit 1322 enddo 1323 enddo 1324 elseif(textpart(1)(1:13).eq.'*FLUIDSECTION') then 1325 nconstants=-1 1326 do i=2,n 1327 if(textpart(i)(1:10).eq.'CONSTANTS=') then 1328 read(textpart(i)(11:20),'(i10)',iostat=istat) 1329 & nconstants 1330 if(istat.gt.0) then 1331 call inputerror(inpc,ipoinpc,iline, 1332 & "*FLUID SECTION%",ier) 1333 exit 1334 endif 1335 nprop_=nprop_+nconstants 1336 exit 1337 endif 1338 enddo 1339 if(nconstants.lt.0) nprop_=nprop_+65 1340 do 1341 call getnewline(inpc,textpart,istat,n,key,iline,ipol,inl, 1342 & ipoinp,inp,ipoinpc) 1343 if((istat.lt.0).or.(key.eq.1)) exit 1344 enddo 1345 elseif(textpart(1)(1:9).eq.'*FRICTION') then 1346! 1347! '8' is for Mortar. 1348! 1349 ncmat_=max(8,ncmat_) 1350 ntmat_=max(1,ntmat_) 1351 call getnewline(inpc,textpart,istat,n,key,iline,ipol,inl, 1352 & ipoinp,inp,ipoinpc) 1353 elseif(textpart(1)(1:5).eq.'*GAP ') then 1354 nmat_=nmat_+1 1355 ncmat_=max(6,ncmat_) 1356 ntmat_=max(1,ntmat_) 1357 do 1358 call getnewline(inpc,textpart,istat,n,key,iline,ipol, 1359 & inl,ipoinp,inp,ipoinpc) 1360 if((istat.lt.0).or.(key.eq.1)) exit 1361 enddo 1362 elseif(textpart(1)(1:15).eq.'*GAPCONDUCTANCE') then 1363 ntmatl=0 1364 do 1365 call getnewline(inpc,textpart,istat,n,key,iline,ipol, 1366 & inl,ipoinp,inp,ipoinpc) 1367 if((istat.lt.0).or.(key.eq.1)) exit 1368 read(textpart(3)(1:20),'(f20.0)',iostat=istat) 1369 & temperature 1370 if(istat.gt.0) then 1371 call inputerror(inpc,ipoinpc,iline, 1372 & "*GAP CONDUCTANCE%",ier) 1373 exit 1374 endif 1375 if(ntmatl.eq.0) then 1376 npmatl=0 1377 ntmatl=ntmatl+1 1378 ntmat_=max(ntmatl,ntmat_) 1379 tempact=temperature 1380 elseif(temperature.ne.tempact) then 1381 npmatl=0 1382 ntmatl=ntmatl+1 1383 ntmat_=max(ntmatl,ntmat_) 1384 tempact=temperature 1385 endif 1386 npmatl=npmatl+1 1387 npmat_=max(npmatl,npmat_) 1388 enddo 1389 elseif(textpart(1)(1:18).eq.'*GAPHEATGENERATION') then 1390 ncmat_=max(11,ncmat_) 1391 call getnewline(inpc,textpart,istat,n,key,iline,ipol,inl, 1392 & ipoinp,inp,ipoinpc) 1393 elseif(textpart(1)(1:20).eq.'*GEOMETRICCONSTRAINT') then 1394 do 1395 call getnewline(inpc,textpart,istat,n,key,iline,ipol, 1396 & inl,ipoinp,inp,ipoinpc) 1397 if((istat.lt.0).or.(key.eq.1)) exit 1398 nobject_ = nobject_ + 1 1399 enddo 1400 elseif(textpart(1)(1:8).eq.'*HEADING') then 1401 if(nheading_.ne.0) then 1402 write(*,*) '*ERROR in allocation: more than 1' 1403 write(*,*) ' *HEADING card in the input deck' 1404 call inputerror(inpc,ipoinpc,iline, 1405 & "*HEADING%",ier) 1406 cycle loop 1407 endif 1408 do 1409 call getnewline(inpc,textpart,istat,n,key,iline,ipol,inl, 1410 & ipoinp,inp,ipoinpc) 1411 if((istat.lt.0).or.(key.eq.1)) exit 1412 nheading_=nheading_+1 1413 enddo 1414 elseif(textpart(1)(1:13).eq.'*HYPERELASTIC') then 1415 ntmatl=0 1416 ityp=-7 1417 do i=2,n 1418 if(textpart(i)(1:12).eq.'ARRUDA-BOYCE') then 1419 ityp=-1 1420 ncmat_=max(3,ncmat_) 1421 elseif(textpart(i)(1:13).eq.'MOONEY-RIVLIN') then 1422 ityp=-2 1423 ncmat_=max(3,ncmat_) 1424 elseif(textpart(i)(1:8).eq.'NEOHOOKE') then 1425 ityp=-3 1426 ncmat_=max(2,ncmat_) 1427 elseif(textpart(i)(1:5).eq.'OGDEN') then 1428 ityp=-4 1429 ncmat_=max(3,ncmat_) 1430 elseif(textpart(i)(1:10).eq.'POLYNOMIAL') then 1431 ityp=-7 1432 ncmat_=max(3,ncmat_) 1433 elseif(textpart(i)(1:17).eq.'REDUCEDPOLYNOMIAL') 1434 & then 1435 ityp=-10 1436 ncmat_=max(2,ncmat_) 1437 elseif(textpart(i)(1:11).eq.'VANDERWAALS') then 1438 ityp=-13 1439 ncmat_=max(5,ncmat_) 1440 elseif(textpart(i)(1:4).eq.'YEOH') then 1441 ityp=-14 1442 ncmat_=max(6,ncmat_) 1443 elseif(textpart(i)(1:2).eq.'N=') then 1444 if(textpart(i)(3:3).eq.'1') then 1445 elseif(textpart(i)(3:3).eq.'2') then 1446 if(ityp.eq.-4) then 1447 ityp=-5 1448 ncmat_=max(6,ncmat_) 1449 elseif(ityp.eq.-7) then 1450 ityp=-8 1451 ncmat_=max(7,ncmat_) 1452 elseif(ityp.eq.-10) then 1453 ityp=-11 1454 ncmat_=max(4,ncmat_) 1455 endif 1456 elseif(textpart(i)(3:3).eq.'3') then 1457 if(ityp.eq.-4) then 1458 ityp=-6 1459 ncmat_=max(9,ncmat_) 1460 elseif(ityp.eq.-7) then 1461 ityp=-9 1462 ncmat_=max(12,ncmat_) 1463 elseif(ityp.eq.-10) then 1464 ityp=-12 1465 ncmat_=max(6,ncmat_) 1466 endif 1467 endif 1468 endif 1469 enddo 1470 if((ityp.ne.-6).and.(ityp.ne.-9)) then 1471 do 1472 call getnewline(inpc,textpart,istat,n,key,iline,ipol, 1473 & inl,ipoinp,inp,ipoinpc) 1474 if((istat.lt.0).or.(key.eq.1)) exit 1475 ntmatl=ntmatl+1 1476 ntmat_=max(ntmatl,ntmat_) 1477 enddo 1478 else 1479 do 1480 call getnewline(inpc,textpart,istat,n,key,iline,ipol, 1481 & inl,ipoinp,inp,ipoinpc) 1482 if((istat.lt.0).or.(key.eq.1)) exit 1483 ntmatl=ntmatl+1 1484 ntmat_=max(ntmatl,ntmat_) 1485 iline=iline+1 1486 enddo 1487 endif 1488 elseif(textpart(1)(1:10).eq.'*HYPERFOAM') then 1489 ntmatl=0 1490 ityp=-15 1491 ncmat_=max(3,ncmat_) 1492 do i=2,n 1493 if(textpart(i)(1:2).eq.'N=') then 1494 if(textpart(i)(3:3).eq.'1') then 1495 elseif(textpart(i)(3:3).eq.'2') then 1496 ityp=-16 1497 ncmat_=max(6,ncmat_) 1498 elseif(textpart(i)(3:3).eq.'3') then 1499 ityp=-17 1500 ncmat_=max(9,ncmat_) 1501 endif 1502 endif 1503 enddo 1504 if(ityp.ne.-17) then 1505 do 1506 call getnewline(inpc,textpart,istat,n,key,iline,ipol, 1507 & inl,ipoinp,inp,ipoinpc) 1508 if((istat.lt.0).or.(key.eq.1)) exit 1509 ntmatl=ntmatl+1 1510 ntmat_=max(ntmatl,ntmat_) 1511 enddo 1512 else 1513 do 1514 call getnewline(inpc,textpart,istat,n,key,iline,ipol, 1515 & inl,ipoinp,inp,ipoinpc) 1516 if((istat.lt.0).or.(key.eq.1)) exit 1517 ntmatl=ntmatl+1 1518 ntmat_=max(ntmatl,ntmat_) 1519 iline=iline+1 1520 enddo 1521 endif 1522 elseif(textpart(1)(2:10).eq.'KINEMATIC') then 1523 npt_=max(npt_,numnodes) 1524! 1525! connection of rotational dofs in refnode to 1526! translational dofs in rotational node 1527! 1528 nk_=nk_+1 1529 nmpc_=nmpc_+3 1530 memmpc_=memmpc_+6 1531! 1532! local system 1533! 1534 if(iorientation.ne.0) then 1535 nk_=nk_+2*numnodes 1536 nmpc_=nmpc_+3*numnodes 1537 memmpc_=memmpc_+3*6*numnodes 1538 nboun_=nboun_+3*numnodes 1539 endif 1540! 1541 do 1542 call getnewline(inpc,textpart,istat,n,key,iline,ipol,inl, 1543 & ipoinp,inp,ipoinpc) 1544 if((istat.lt.0).or.(key.eq.1)) exit 1545! 1546 read(textpart(1)(1:10),'(i10)',iostat=istat) ibounstart 1547 if(istat.gt.0) then 1548 call inputerror(inpc,ipoinpc,iline, 1549 & "*BOUNDARY%",ier) 1550 exit 1551 endif 1552! 1553 if(textpart(2)(1:1).eq.' ') then 1554 ibounend=ibounstart 1555 else 1556 read(textpart(2)(1:10),'(i10)',iostat=istat) ibounend 1557 if(istat.gt.0) then 1558 call inputerror(inpc,ipoinpc,iline, 1559 & "*BOUNDARY%",ier) 1560 exit 1561 endif 1562 endif 1563 ibound=ibounend-ibounstart+1 1564 ibound=max(1,ibound) 1565 ibound=min(3,ibound) 1566! 1567 if(iorientation.eq.0) then 1568 nk_=nk_+numnodes 1569 nmpc_=nmpc_+ibound*numnodes 1570 memmpc_=memmpc_+6*ibound*numnodes 1571 nboun_=nboun_+ibound*numnodes 1572 else 1573 nmpc_=nmpc_+ibound*numnodes 1574 memmpc_=memmpc_+ibound*6*numnodes 1575 endif 1576 enddo 1577 elseif(textpart(1)(1:21).eq.'*MAGNETICPERMEABILITY') then 1578 ntmatl=0 1579 ityp=2 1580 ncmat_=max(2,ncmat_) 1581 do 1582 call getnewline(inpc,textpart,istat,n,key,iline,ipol, 1583 & inl,ipoinp,inp,ipoinpc) 1584 if((istat.lt.0).or.(key.eq.1)) exit 1585 ntmatl=ntmatl+1 1586 ntmat_=max(ntmatl,ntmat_) 1587 enddo 1588 elseif(textpart(1)(1:5).eq.'*MASS') then 1589 nmat_=nmat_+1 1590 ntmat_=max(1,ntmat_) 1591 call getnewline(inpc,textpart,istat,n,key,iline,ipol, 1592 & inl,ipoinp,inp,ipoinpc) 1593 elseif(textpart(1)(1:9).eq.'*MATERIAL') then 1594 nmat_=nmat_+1 1595 call getnewline(inpc,textpart,istat,n,key,iline,ipol,inl, 1596 & ipoinp,inp,ipoinpc) 1597 elseif(textpart(1)(1:13).eq.'*MODALDAMPING') then 1598 if(textpart(2)(1:8).ne.'RAYLEIGH') then 1599 nevdamp_=0 1600 do 1601 call getnewline(inpc,textpart,istat,n,key,iline,ipol, 1602 & inl,ipoinp,inp,ipoinpc) 1603 if((istat.lt.0).or.(key.eq.1)) exit 1604 read(textpart(1)(1:10),'(i10)',iostat=istat) i 1605 if(istat.gt.0) then 1606 call inputerror(inpc,ipoinpc,iline, 1607 & "*MODAL DAMPING%",ier) 1608 exit 1609 endif 1610 nevdamp_ = max(nevdamp_,i) 1611 read(textpart(2)(1:10),'(i10)',iostat=istat) i 1612 if(istat.gt.0) then 1613 call inputerror(inpc,ipoinpc,iline, 1614 & "*MODAL DAMPING%",ier) 1615 exit 1616 endif 1617 nevdamp_ = max(nevdamp_,i) 1618 enddo 1619 else 1620 call getnewline(inpc,textpart,istat,n,key,iline,ipol,inl, 1621 & ipoinp,inp,ipoinpc) 1622 endif 1623 elseif(textpart(1)(1:12).eq.'*MODELCHANGE') then 1624 if(iprestr.ne.2) then 1625 do i=2,n 1626 if(textpart(i)(1:14).eq.'ADD=STRAINFREE') then 1627 iprestr=2 1628 exit 1629 elseif(textpart(i)(1:14).eq.'ADD=WITHSTRAIN') then 1630 elseif(textpart(i)(1:3).eq.'ADD') then 1631 iprestr=2 1632 exit 1633 elseif(textpart(i)(1:20).eq.'MECHSTRAINTORESIDUAL') 1634 & then 1635 iprestr=2 1636 exit 1637 endif 1638 enddo 1639 endif 1640 call getnewline(inpc,textpart,istat,n,key,iline,ipol,inl, 1641 & ipoinp,inp,ipoinpc) 1642 elseif(textpart(1)(1:4).eq.'*MPC') then 1643 mpclabel=' ' 1644 do 1645 call getnewline(inpc,textpart,istat,n,key,iline,ipol,inl, 1646 & ipoinp,inp,ipoinpc) 1647 if((istat.lt.0).or.(key.eq.1)) exit 1648 do i=1,n 1649 read(textpart(i)(1:10),'(i10)',iostat=istat) ialset(i) 1650 if(mpclabel.eq.' ') then 1651 mpclabel=textpart(i)(1:20) 1652 if((mpclabel(1:8).ne.'STRAIGHT').and. 1653 & (mpclabel(1:4).ne.'PLANE')) then 1654 nk_=nk_+1 1655 nmpc_=nmpc_+1 1656 nboun_=nboun_+1 1657 memmpc_=memmpc_+1 1658 endif 1659 elseif(istat.gt.0) then 1660 noelset=textpart(i)(1:80) 1661 noelset(81:81)=' ' 1662 ipos=index(noelset,' ') 1663 noelset(ipos:ipos)='N' 1664c do j=1,nset_ 1665c if(noelset.eq.set(j)) then 1666 call cident81(set,noelset,nset_,j) 1667 if(j.gt.0) then 1668 if(set(j).eq.noelset) then 1669 if(mpclabel(1:8).eq.'STRAIGHT') then 1670 nk_=nk_+2*meminset(j) 1671 nmpc_=nmpc_+2*meminset(j) 1672 nboun_=nboun_+2*meminset(j) 1673 memmpc_=memmpc_+14*meminset(j) 1674 elseif(mpclabel(1:5).eq.'PLANE') then 1675 nk_=nk_+meminset(j) 1676 nmpc_=nmpc_+meminset(j) 1677 nboun_=nboun_+meminset(j) 1678 memmpc_=memmpc_+13*meminset(j) 1679 elseif(mpclabel(1:4).eq.'BEAM') then 1680 memmpc_=memmpc_+3*meminset(j) 1681 else 1682 memmpc_=memmpc_+meminset(j) 1683 endif 1684c exit 1685 endif 1686 endif 1687 else 1688 if(mpclabel(1:8).eq.'STRAIGHT') then 1689 nk_=nk_+2 1690 nmpc_=nmpc_+2 1691 nboun_=nboun_+2 1692 memmpc_=memmpc_+14 1693 elseif(mpclabel(1:5).eq.'PLANE') then 1694 nk_=nk_+1 1695 nmpc_=nmpc_+1 1696 nboun_=nboun_+1 1697 memmpc_=memmpc_+13 1698 elseif(mpclabel(1:4).eq.'BEAM') then 1699 memmpc_=memmpc_+3 1700 else 1701 memmpc_=memmpc_+1 1702 endif 1703 endif 1704 enddo 1705 enddo 1706 elseif(textpart(1)(1:11).eq.'*NETWORKMPC') then 1707 do 1708 call getnewline(inpc,textpart,istat,n,key,iline,ipol,inl, 1709 & ipoinp,inp,ipoinpc) 1710 if((istat.lt.0).or.(key.eq.1)) exit 1711 read(textpart(1)(1:10),'(i10)',iostat=istat) nterm 1712 if(ntrans_.eq.0) then 1713 nmpc_=nmpc_+1 1714 memmpc_=memmpc_+nterm 1715 else 1716 nmpc_=nmpc_+3 1717 memmpc_=memmpc_+3*nterm 1718 endif 1719 ii=0 1720 do 1721 call getnewline(inpc,textpart,istat,n,key,iline,ipol, 1722 & inl,ipoinp,inp,ipoinpc) 1723 if((istat.lt.0).or.(key.eq.1)) exit 1724 ii=ii+n/3 1725 if(ii.eq.nterm) exit 1726 enddo 1727 enddo 1728 elseif((textpart(1)(1:5).eq.'*NODE').and. 1729 & (textpart(1)(1:10).ne.'*NODEPRINT').and. 1730 & (textpart(1)(1:9).ne.'*NODEFILE').and. 1731 & (textpart(1)(1:11).ne.'*NODEOUTPUT')) then 1732 inoset=0 1733 loop3: do i=2,n 1734 if(textpart(i)(1:5).eq.'NSET=') then 1735 noset=textpart(i)(6:85) 1736 noset(81:81)=' ' 1737 ipos=index(noset,' ') 1738 noset(ipos:ipos)='N' 1739 inoset=1 1740c do js=1,nset_ 1741c if(set(js).eq.noset) exit 1742c enddo 1743c if(js.gt.nset_) then 1744c nset_=nset_+1 1745c set(nset_)=noset 1746c endif 1747 call cident81(set,noset,nset_,id) 1748 js=nset_+1 1749 if(id.gt.0) then 1750 if(set(id).eq.noset) js=id 1751 endif 1752 if(js.gt.nset_) then 1753 nset_=nset_+1 1754 do j=nset_,id+2,-1 1755 meminset(j)=meminset(j-1) 1756 rmeminset(j)=rmeminset(j-1) 1757 set(j)=set(j-1) 1758 enddo 1759 js=id+1 1760 set(js)=noset 1761 meminset(js)=0 1762 rmeminset(js)=0 1763 endif 1764 endif 1765 enddo loop3 1766! 1767 do 1768 call getnewline(inpc,textpart,istat,n,key,iline,ipol,inl, 1769 & ipoinp,inp,ipoinpc) 1770 if((istat.lt.0).or.(key.eq.1)) exit 1771 read(textpart(1)(1:10),'(i10)',iostat=istat) i 1772 if(istat.gt.0) then 1773 call inputerror(inpc,ipoinpc,iline, 1774 & "*NODE or *NODE PRINT or *NODE FILE or *NODE OUTPUT%", 1775 & ier) 1776 exit 1777 endif 1778 nk_=max(nk_,i) 1779 if(inoset.eq.1) then 1780 meminset(js)=meminset(js)+1 1781 rmeminset(js)=rmeminset(js)+1 1782 endif 1783 enddo 1784 elseif(textpart(1)(1:10).eq.'*NODEPRINT') then 1785 do 1786 call getnewline(inpc,textpart,istat,n,key,iline,ipol,inl, 1787 & ipoinp,inp,ipoinpc) 1788 if((istat.lt.0).or.(key.eq.1)) exit 1789 nprint_=nprint_+n 1790 enddo 1791 elseif(textpart(1)(1:10).eq.'*OBJECTIVE') then 1792 nobject_=nobject_+1 1793 call getnewline(inpc,textpart,istat,n,key,iline,ipol,inl, 1794 & ipoinp,inp,ipoinpc) 1795 elseif(textpart(1)(1:12).eq.'*ORIENTATION') then 1796 norien_=norien_+1 1797 do 1798 call getnewline(inpc,textpart,istat,n,key,iline,ipol,inl, 1799 & ipoinp,inp,ipoinpc) 1800 if((istat.lt.0).or.(key.eq.1)) exit 1801 enddo 1802 elseif(textpart(1)(1:8).eq.'*PLASTIC') then 1803 ntmatl=0 1804 do 1805 call getnewline(inpc,textpart,istat,n,key,iline,ipol,inl, 1806 & ipoinp,inp,ipoinpc) 1807 if((istat.lt.0).or.(key.eq.1)) exit 1808 read(textpart(3)(1:20),'(f20.0)',iostat=istat) 1809 & temperature 1810 if(istat.gt.0) then 1811 call inputerror(inpc,ipoinpc,iline, 1812 & "*PLASTIC%",ier) 1813 exit 1814 endif 1815 if(ntmatl.eq.0) then 1816 npmatl=0 1817 ntmatl=ntmatl+1 1818 ntmat_=max(ntmatl,ntmat_) 1819 tempact=temperature 1820 elseif(temperature.ne.tempact) then 1821 npmatl=0 1822 ntmatl=ntmatl+1 1823 ntmat_=max(ntmatl,ntmat_) 1824 tempact=temperature 1825 endif 1826 npmatl=npmatl+1 1827 npmat_=max(npmatl,npmat_) 1828 enddo 1829 if(ncmat_.ge.9) ncmat_=max(19,ncmat_) 1830 elseif(textpart(1)(1:19).eq.'*PRE-TENSIONSECTION') then 1831 surface(1:1)=' ' 1832 do i=2,n 1833 if(textpart(i)(1:8).eq.'SURFACE=') then 1834 surface=textpart(i)(9:88) 1835 ipos=index(surface,' ') 1836 surface(ipos:ipos)='T' 1837 exit 1838 elseif(textpart(i)(1:8).eq.'ELEMENT=') then 1839 nmpc_=nmpc_+1 1840 memmpc_=memmpc_+7 1841 exit 1842 endif 1843 enddo 1844 if(surface(1:1).ne.' ') then 1845c do i=1,nset_ 1846c if(set(i).eq.surface) then 1847 call cident81(set,surface,nset_,i) 1848 if(i.gt.0) then 1849 if(set(i).eq.surface) then 1850! 1851! worst case: 8 nodes per element face 1852! 1853 nk_=nk_+8*meminset(i) 1854 npt_=npt_+8*meminset(i) 1855! 1856! 2 MPC's per node perpendicular to tension direction 1857! + 1 thermal MPC per node 1858! + 1 MPC per node in tension direction (the total of 1859! which is divided into one global tension MPC and the 1860! rest are MPC's specifying that the distance in tension 1861! direction in all nodes should be the same) 1862! 1863 nmpc_=nmpc_+32*meminset(i)+1 1864! 1865! 6 terms per MPC perpendicular to tension direction 1866! + 2 thermal terms per MPC 1867! + 6 terms * # of nodes +1 parallel to tension 1868! direction 1869! + 12 terms per MPC parallel to tension direction 1870! 1871 memmpc_=memmpc_+96*meminset(i) 1872 & +16*meminset(i) 1873 & +48*meminset(i)+1 1874 & +12*(8*meminset(i)-1) 1875c exit 1876! 1877 endif 1878 endif 1879 endif 1880 call getnewline(inpc,textpart,istat,n,key,iline,ipol,inl, 1881 & ipoinp,inp,ipoinpc) 1882 elseif(textpart(1)(1:8).eq.'*RADIATE') then 1883 nam_=nam_+2 1884 namtot_=namtot_+2 1885 do 1886 call getnewline(inpc,textpart,istat,n,key,iline,ipol,inl, 1887 & ipoinp,inp,ipoinpc) 1888 if((istat.lt.0).or.(key.eq.1)) exit 1889 read(textpart(2)(1:5),'(a5)',iostat=istat) llab 1890 if((llab.eq.'GRAV ').or.(llab.eq.'CENTR')) exit 1891 read(textpart(1)(1:10),'(i10)',iostat=istat) l 1892 if(istat.eq.0) then 1893 nload_=nload_+1 1894 nradiate=nradiate+1 1895 else 1896 read(textpart(1)(1:80),'(a80)',iostat=istat) elset 1897 elset(81:81)=' ' 1898 ipos=index(elset,' ') 1899 elset(ipos:ipos)='E' 1900c do i=1,nset_ 1901c if(set(i).eq.elset) then 1902 call cident81(set,elset,nset_,i) 1903 if(i.gt.0) then 1904 if(set(i).eq.elset) then 1905 nload_=nload_+meminset(i) 1906 nradiate=nradiate+meminset(i) 1907c exit 1908 endif 1909 endif 1910 endif 1911 enddo 1912 elseif(textpart(1)(1:8).eq.'*RESTART') then 1913 irestartread=0 1914 irestartstep=0 1915 do i=1,n 1916 if(textpart(i)(1:4).eq.'READ') then 1917 irestartread=1 1918 endif 1919 if(textpart(i)(1:5).eq.'STEP=') then 1920 read(textpart(i)(6:15),'(i10)',iostat=istat) 1921 & irestartstep 1922 if(istat.gt.0) then 1923 call inputerror(inpc,ipoinpc,iline, 1924 & "*RESTART%",ier) 1925 exit 1926 endif 1927 endif 1928 enddo 1929 if(irestartread.eq.1) then 1930 icntrl=1 1931 call restartshort(nset_,nload_,nbody_,nforc_,nboun_,nk_,ne_, 1932 & nmpc_,nalset_,nmat_,ntmat_,npmat_,norien_,nam_,nprint_, 1933 & mi,ntrans_,ncs_,namtot_,ncmat_,memmpc_, 1934 & ne1d,ne2d,nflow,set,meminset,rmeminset,jobnamec, 1935 & irestartstep,icntrl,ithermal,nener,nstate_,ntie_, 1936 & nslavs,nkon_,mcs,nprop_,mortar,ifacecount,nintpoint, 1937 & infree,nef,mpcend,nheading_,network) 1938 irstrt(1)=-1 1939 nbounold=nboun_ 1940 nforcold=nforc_ 1941 nloadold=nload_ 1942 nbodyold=nbody_ 1943 else 1944 endif 1945 call getnewline(inpc,textpart,istat,n,key,iline,ipol,inl, 1946 & ipoinp,inp,ipoinpc) 1947 elseif(textpart(1)(1:18).eq.'*RETAINEDNODALDOFS') then 1948 do 1949 call getnewline(inpc,textpart,istat,n,key,iline,ipol,inl, 1950 & ipoinp,inp,ipoinpc) 1951 if((istat.lt.0).or.(key.eq.1)) exit 1952! 1953 read(textpart(2)(1:10),'(i10)',iostat=istat) ibounstart 1954 if(istat.gt.0) then 1955 call inputerror(inpc,ipoinpc,iline, 1956 & "*BOUNDARY%",ier) 1957 exit 1958 endif 1959! 1960 if(textpart(3)(1:1).eq.' ') then 1961 ibounend=ibounstart 1962 else 1963 read(textpart(3)(1:10),'(i10)',iostat=istat) ibounend 1964 if(istat.gt.0) then 1965 call inputerror(inpc,ipoinpc,iline, 1966 & "*BOUNDARY%",ier) 1967 exit 1968 endif 1969 endif 1970 ibound=ibounend-ibounstart+1 1971 ibound=max(1,ibound) 1972 ibound=min(3,ibound) 1973! 1974 read(textpart(1)(1:10),'(i10)',iostat=istat) l 1975 if(istat.eq.0) then 1976 nboun_=nboun_+ibound 1977 else 1978 read(textpart(1)(1:80),'(a80)',iostat=istat) noset 1979 noset(81:81)=' ' 1980 ipos=index(noset,' ') 1981 noset(ipos:ipos)='N' 1982c do i=1,nset_ 1983c if(set(i).eq.noset) then 1984 call cident81(set,noset,nset_,i) 1985 if(i.gt.0) then 1986 if(set(i).eq.noset) then 1987 nboun_=nboun_+ibound*meminset(i) 1988c exit 1989 endif 1990 endif 1991 endif 1992 enddo 1993 elseif(textpart(1)(1:10).eq.'*RIGIDBODY') then 1994 noset=' 1995 &' 1996 elset=' 1997 &' 1998 do i=2,n 1999 if(textpart(i)(1:5).eq.'NSET=') 2000 & then 2001 noset=textpart(i)(6:85) 2002 noset(81:81)=' ' 2003 ipos=index(noset,' ') 2004 noset(ipos:ipos)='N' 2005 exit 2006 elseif(textpart(i)(1:6).eq.'ELSET=') 2007 & then 2008 elset=textpart(i)(7:86) 2009 elset(81:81)=' ' 2010 ipos=index(elset,' ') 2011 elset(ipos:ipos)='E' 2012 exit 2013 endif 2014 enddo 2015 if(noset(1:1).ne.' ') then 2016c do i=1,nset_ 2017c if(set(i).eq.noset) then 2018 call cident81(set,noset,nset_,i) 2019 if(i.gt.0) then 2020 if(set(i).eq.noset) then 2021 nk_=nk_+2+meminset(i) 2022 nmpc_=nmpc_+3*meminset(i) 2023 memmpc_=memmpc_+18*meminset(i) 2024 nboun_=nboun_+3*meminset(i) 2025 endif 2026 endif 2027 elseif(elset(1:1).ne.' ') then 2028c do i=1,nset_ 2029c if(set(i).eq.elset) then 2030 call cident81(set,elset,nset_,i) 2031 if(i.gt.0) then 2032 if(set(i).eq.elset) then 2033 nk_=nk_+2+20*meminset(i) 2034 nmpc_=nmpc_+60*meminset(i) 2035 memmpc_=memmpc_+360*meminset(i) 2036 nboun_=nboun_+60*meminset(i) 2037 endif 2038 endif 2039 endif 2040 call getnewline(inpc,textpart,istat,n,key,iline,ipol,inl, 2041 & ipoinp,inp,ipoinpc) 2042 elseif(textpart(1)(1:13).eq.'*ROBUSTDESIGN') then 2043 call getnewline(inpc,textpart,istat,n,key,iline,ipol,inl, 2044 & ipoinp,inp,ipoinpc) 2045 if((istat.lt.0).or.(key.eq.1)) exit 2046 irobustdesign(1)=1 2047 elseif(textpart(1)(1:16).eq.'*SECTIONPRINT') then 2048 do 2049 call getnewline(inpc,textpart,istat,n,key,iline,ipol,inl, 2050 & ipoinp,inp,ipoinpc) 2051 if((istat.lt.0).or.(key.eq.1)) exit 2052 nprint_=nprint_+n 2053 enddo 2054 elseif(textpart(1)(1:13).eq.'*SHELLSECTION') then 2055 composite=.false. 2056 do i=2,n 2057 if(textpart(i)(1:9).eq.'COMPOSITE') then 2058 composite=.true. 2059 nlayer=0 2060 elseif(textpart(i)(1:6).eq.'ELSET=') then 2061 elset=textpart(i)(7:86) 2062 elset(81:81)=' ' 2063 ipos=index(elset,' ') 2064 elset(ipos:ipos)='E' 2065c do js=1,nset_ 2066c if(set(js).eq.elset) exit 2067c enddo 2068 call cident81(set,elset,nset_,id) 2069 js=nset_+1 2070 if(id.gt.0) then 2071 if(set(id).eq.elset) then 2072 js=id 2073 endif 2074 endif 2075 endif 2076 enddo 2077 if(composite) then 2078 do 2079 call getnewline(inpc,textpart,istat,n,key,iline,ipol, 2080 & inl,ipoinp,inp,ipoinpc) 2081 if((istat.lt.0).or.(key.eq.1)) then 2082! 2083! conservative upper limit 2084! "label" is not necessary the label of the 2085! composite shell element 2086! 2087 mi(1)=max(mi(1),8*nlayer) 2088 mi(3)=max(mi(3),nlayer) 2089 if(js.le.nset_) then 2090 nk_=nk_+20*nlayer*meminset(js) 2091 nkon_=nkon_+20*nlayer*meminset(js) 2092 endif 2093 exit 2094 endif 2095 nlayer=nlayer+1 2096 enddo 2097 else 2098 call getnewline(inpc,textpart,istat,n,key,iline,ipol, 2099 & inl,ipoinp,inp,ipoinpc) 2100 endif 2101 elseif(textpart(1)(1:7).eq.'*SPRING') then 2102 nmat_=nmat_+1 2103 lin=.true. 2104 do i=2,n 2105 if(textpart(i)(1:9).eq.'NONLINEAR') then 2106 lin=.false. 2107 exit 2108 endif 2109 enddo 2110 if(lin) then 2111 ntmatl=0 2112 ncmat_=max(2,ncmat_) 2113 do 2114 call getnewline(inpc,textpart,istat,n,key,iline,ipol, 2115 & inl,ipoinp,inp,ipoinpc) 2116 if((istat.lt.0).or.(key.eq.1)) exit 2117 ntmatl=ntmatl+1 2118 ntmat_=max(ntmatl,ntmat_) 2119 enddo 2120 else 2121 ntmatl=0 2122 do 2123 call getnewline(inpc,textpart,istat,n,key,iline,ipol, 2124 & inl,ipoinp,inp,ipoinpc) 2125 if((istat.lt.0).or.(key.eq.1)) exit 2126 read(textpart(3)(1:20),'(f20.0)',iostat=istat) 2127 & temperature 2128 if(istat.gt.0) then 2129 call inputerror(inpc,ipoinpc,iline, 2130 & "*SPRING%",ier) 2131 exit 2132 endif 2133 if(ntmatl.eq.0) then 2134 npmatl=0 2135 ntmatl=ntmatl+1 2136 ntmat_=max(ntmatl,ntmat_) 2137 tempact=temperature 2138 elseif(temperature.ne.tempact) then 2139 npmatl=0 2140 ntmatl=ntmatl+1 2141 ntmat_=max(ntmatl,ntmat_) 2142 tempact=temperature 2143 endif 2144 npmatl=npmatl+1 2145 npmat_=max(npmatl,npmat_) 2146 enddo 2147 if(ncmat_.ge.9) ncmat_=max(19,ncmat_) 2148 endif 2149 elseif(textpart(1)(1:5).eq.'*STEP') then 2150 if(nstam.eq.0) then 2151 do i=1,n 2152 if((textpart(i)(1:14).eq.'AMPLITUDE=STEP').or. 2153 & (textpart(i)(1:14).eq.'AMPLITUDE=RAMP')) then 2154 nam_=nam_+2 2155 namtot_=namtot_+4 2156 nstam=1 2157 exit 2158 endif 2159 enddo 2160 endif 2161 call getnewline(inpc,textpart,istat,n,key,iline,ipol,inl, 2162 & ipoinp,inp,ipoinpc) 2163 elseif(textpart(1)(1:9).eq.'*SUBMODEL') then 2164 nsubmodel=nsubmodel+1 2165 ntie_=ntie_+1 2166 nam_=nam_+1 2167 namtot_=namtot_+4 2168! 2169! global element set 2170! 2171 do j=2,n 2172 if(textpart(j)(1:12).eq.'GLOBALELSET=') 2173 & then 2174 mastset(1:80)=textpart(j)(13:92) 2175 mastset(81:81)=' ' 2176 ipos=index(mastset,' ') 2177 mastset(ipos:ipos)='E' 2178c do i=1,nset_ 2179c if(set(i).eq.mastset) exit 2180c enddo 2181 call cident81(set,mastset,nset_,id) 2182 i=nset_+1 2183 if(id.gt.0) then 2184 if(set(id).eq.mastset) then 2185 i=id 2186 endif 2187 endif 2188 if(i.le.nset_) then 2189cc nset_=nset_+1 2190cc do k=1,81 2191cc set(nset_)(k:k)=' ' 2192cc enddo 2193c globalset(1:7)=' GLOBAL' 2194c do k=8,81 2195c globalset(k:k)=' ' 2196c enddo 2197c memglobal=meminset(i) 2198c call cident81(set,globalset,nset_,id) 2199c nset_=nset_+1 2200c do k=nset_,id+2,-1 2201c meminset(k)=meminset(k-1) 2202c rmeminset(k)=rmeminset(k-1) 2203c set(k)=set(k-1) 2204c enddo 2205c set(id+1)=globalset 2206c meminset(id+1)=meminset(id+1)+memglobal 2207c rmeminset(id+1)=rmeminset(id+1)+memglobal 2208cc meminset(nset_)=meminset(nset_)+meminset(i) 2209cc rmeminset(nset_)=rmeminset(nset_)+meminset(i) 2210 meminset(i)=meminset(i)+meminset(i) 2211 rmeminset(i)=rmeminset(i)+meminset(i) 2212 endif 2213 elseif(textpart(j)(1:5).eq.'TYPE=') then 2214 if(textpart(j)(6:12).eq.'SURFACE') then 2215 selabel='T' 2216 else 2217 selabel='N' 2218 endif 2219 endif 2220 enddo 2221! 2222! local node or element face set 2223! 2224c nset_=nset_+1 2225c submset(1:6)=' LOCAL' 2226c do k=7,81 2227c globalset(k:k)=' ' 2228c enddo 2229 submset(1:8)='SUBMODEL' 2230 if(nsubmodel.lt.10) then 2231 submset(9:10)='00' 2232 write(submset(11:11),'(i1)') nsubmodel 2233 elseif(nsubmodel.lt.100) then 2234 submset(9:9)='0' 2235 write(submset(10:11),'(i2)') nsubmodel 2236 elseif(nsubmodel.lt.1000) then 2237 write(submset(9:11),'(i3)') nsubmodel 2238 else 2239 write(*,*) '*ERROR reading *SUBMODEL: no more than 999' 2240 write(*,*) ' submodels allowed' 2241 ier=1 2242 return 2243 endif 2244 submset(12:12)=selabel 2245 do i=13,81 2246 submset(i:i)=' ' 2247 enddo 2248! 2249 call cident81(set,submset,nset_,id) 2250 nset_=nset_+1 2251 do j=nset_,id+2,-1 2252 meminset(j)=meminset(j-1) 2253 rmeminset(j)=rmeminset(j-1) 2254 set(j)=set(j-1) 2255 enddo 2256 js=id+1 2257 set(js)=submset 2258 meminset(js)=0 2259 rmeminset(js)=0 2260 do 2261 call getnewline(inpc,textpart,istat,n,key,iline,ipol,inl, 2262 & ipoinp,inp,ipoinpc) 2263 if((istat.lt.0).or.(key.eq.1)) exit 2264 read(textpart(1)(1:10),'(i10)',iostat=istat) ialset(1) 2265 if(istat.gt.0) then 2266 noset=textpart(1)(1:80) 2267 noset(81:81)=' ' 2268 ipos=index(noset,' ') 2269 noset(ipos:ipos)=selabel 2270c do i=1,nset_-1 2271c do i=1,nset_ 2272c if(set(i).eq.noset) then 2273 call cident81(set,noset,nset_,i) 2274 if(i.gt.0) then 2275 if(set(i).eq.noset) then 2276c meminset(nset_)=meminset(nset_)+meminset(i) 2277 meminset(js)=meminset(js)+meminset(i) 2278! 2279! surfaces are stored in expanded form 2280! (no equivalent to generate) 2281! 2282c rmeminset(nset_)=rmeminset(nset_)+meminset(i) 2283 rmeminset(js)=rmeminset(js)+meminset(i) 2284 endif 2285 endif 2286 else 2287 meminset(js)=meminset(js)+1 2288 rmeminset(js)=rmeminset(js)+1 2289 endif 2290 enddo 2291 elseif(textpart(1)(1:9).eq.'*SURFACE ') then 2292c nset_=nset_+1 2293 sulabel='T' 2294 do i=2,n 2295 if(textpart(i)(1:5).eq.'NAME=') 2296 & then 2297c set(nset_)=textpart(i)(6:85) 2298c set(nset_)(81:81)=' ' 2299 surfset=textpart(i)(6:85) 2300 surfset(81:81)=' ' 2301 elseif(textpart(i)(1:9).eq.'TYPE=NODE') then 2302 sulabel='S' 2303 endif 2304 enddo 2305c ipos=index(set(nset_),' ') 2306c set(nset_)(ipos:ipos)=sulabel 2307 ipos=index(surfset,' ') 2308 surfset(ipos:ipos)=sulabel 2309 call cident81(set,surfset,nset_,id) 2310 nset_=nset_+1 2311 do j=nset_,id+2,-1 2312 meminset(j)=meminset(j-1) 2313 rmeminset(j)=rmeminset(j-1) 2314 set(j)=set(j-1) 2315 enddo 2316 js=id+1 2317 set(js)=surfset 2318 meminset(js)=0 2319 rmeminset(js)=0 2320! 2321 if(sulabel.eq.'S') then 2322 selabel='N' 2323 else 2324 selabel='E' 2325 endif 2326 do 2327 call getnewline(inpc,textpart,istat,n,key,iline,ipol,inl, 2328 & ipoinp,inp,ipoinpc) 2329 if((istat.lt.0).or.(key.eq.1)) exit 2330 read(textpart(1)(1:10),'(i10)',iostat=istat) ialset(1) 2331 if(istat.gt.0) then 2332 noset=textpart(1)(1:80) 2333 noset(81:81)=' ' 2334 ipos=index(noset,' ') 2335 noset(ipos:ipos)=selabel 2336c do i=1,nset_-1 2337c if(set(i).eq.noset) then 2338 call cident81(set,noset,nset_,i) 2339 if(i.gt.0) then 2340 if(set(i).eq.noset) then 2341 meminset(js)=meminset(js)+meminset(i) 2342! 2343! surfaces are stored in expanded form 2344! (no equivalent to generate) 2345! 2346 rmeminset(js)=rmeminset(js)+meminset(i) 2347 endif 2348 endif 2349 else 2350 meminset(js)=meminset(js)+1 2351 rmeminset(js)=rmeminset(js)+1 2352 endif 2353 enddo 2354! 2355! for CFD-calculations: local coordinate systems are 2356! stored as distributed load 2357! 2358 if(icfd>0) nload_=nload_+rmeminset(js) 2359 elseif(textpart(1)(1:16).eq.'*SURFACEBEHAVIOR') then 2360 ncmat_=max(4,ncmat_) 2361 ntmat_=max(1,ntmat_) 2362 tabular=.false. 2363 do i=1,n 2364 if(textpart(i)(1:38).eq.'PRESSURE-OVERCLOSURE=TABULAR') 2365 & tabular=.true. 2366 enddo 2367 if(tabular) then 2368 ntmatl=0 2369 do 2370 call getnewline(inpc,textpart,istat,n,key,iline, 2371 & ipol,inl,ipoinp,inp,ipoinpc) 2372 if((istat.lt.0).or.(key.eq.1)) exit 2373 read(textpart(3)(1:20),'(f20.0)',iostat=istat) 2374 & temperature 2375 if(istat.gt.0) then 2376 call inputerror(inpc,ipoinpc,iline, 2377 & "*SURFACE BEHAVIOR%",ier) 2378 exit 2379 endif 2380 if(ntmatl.eq.0) then 2381 npmatl=0 2382 ntmatl=ntmatl+1 2383 ntmat_=max(ntmatl,ntmat_) 2384 tempact=temperature 2385 elseif(temperature.ne.tempact) then 2386 npmatl=0 2387 ntmatl=ntmatl+1 2388 ntmat_=max(ntmatl,ntmat_) 2389 tempact=temperature 2390 endif 2391 npmatl=npmatl+1 2392 npmat_=max(npmatl,npmat_) 2393 enddo 2394 else 2395 call getnewline(inpc,textpart,istat,n,key,iline,ipol,inl, 2396 & ipoinp,inp,ipoinpc) 2397 endif 2398 elseif(textpart(1)(1:19).eq.'*SURFACEINTERACTION') then 2399 nmat_=nmat_+1 2400 call getnewline(inpc,textpart,istat,n,key,iline,ipol,inl, 2401 & ipoinp,inp,ipoinpc) 2402 elseif(textpart(1)(1:12).eq.'*TEMPERATURE') then 2403 nam_=nam_+1 2404 namtot_=namtot_+1 2405 call getnewline(inpc,textpart,istat,n,key,iline,ipol,inl, 2406 & ipoinp,inp,ipoinpc) 2407 elseif(textpart(1)(1:4).eq.'*TIE') then 2408 ntie_=ntie_+1 2409 cyclicsymmetry=.false. 2410 do i=1,n 2411 if((textpart(i)(1:14).eq.'CYCLICSYMMETRY').or. 2412 & (textpart(i)(1:10).eq.'MULTISTAGE')) then 2413 cyclicsymmetry=.true. 2414 endif 2415 enddo 2416 call getnewline(inpc,textpart,istat,n,key,iline,ipol,inl, 2417 & ipoinp,inp,ipoinpc) 2418 if(.not.cyclicsymmetry) cycle 2419 if((istat.lt.0).or.(key.eq.1)) cycle 2420! 2421 slavset=textpart(1)(1:80) 2422 slavset(81:81)=' ' 2423 iposs=index(slavset,' ') 2424 slavsets=slavset 2425 slavsets(iposs:iposs)='S' 2426 slavsett=slavset 2427 slavsett(iposs:iposs)='T' 2428! 2429 mastset=textpart(2)(1:80) 2430 mastset(81:81)=' ' 2431 iposm=index(mastset,' ') 2432 mastsets=mastset 2433 mastsets(iposm:iposm)='S' 2434 mastsett=mastset 2435 mastsett(iposm:iposm)='T' 2436! 2437 islavset=0 2438 imastset=0 2439! 2440c do i=1,nset_ 2441c if(set(i).eq.slavsets) then 2442c islavset=i 2443c multslav=1 2444c elseif(set(i).eq.slavsett) then 2445c islavset=i 2446c multslav=8 2447c elseif(set(i).eq.mastsets) then 2448c imastset=i 2449c multmast=1 2450c elseif(set(i).eq.mastsett) then 2451c imastset=i 2452c multmast=8 2453c endif 2454c enddo 2455 call cident81(set,slavsets,nset_,i) 2456 if(i.gt.0) then 2457 if(set(i).eq.slavsets) then 2458 islavset=i 2459 multslav=1 2460 endif 2461 endif 2462 call cident81(set,slavsett,nset_,i) 2463 if(i.gt.0) then 2464 if(set(i).eq.slavsett) then 2465 islavset=i 2466 multslav=8 2467 endif 2468 endif 2469 call cident81(set,mastsets,nset_,i) 2470 if(i.gt.0) then 2471 if(set(i).eq.mastsets) then 2472 imastset=i 2473 multmast=1 2474 endif 2475 endif 2476 call cident81(set,mastsett,nset_,i) 2477 if(i.gt.0) then 2478 if(set(i).eq.mastsett) then 2479 imastset=i 2480 multmast=8 2481 endif 2482 endif 2483! 2484 if((islavset.ne.0).and.(imastset.ne.0)) then 2485 ncs_=ncs_+max(multslav*meminset(islavset), 2486 & multmast*meminset(imastset)) 2487 else 2488 write(*,*) '*ERROR in allocation: either the slave' 2489 write(*,*) ' surface or the master surface in a' 2490 write(*,*) ' cyclic symmetry *TIE option or both' 2491 write(*,*) ' do not exist or are no nodal surfaces' 2492 write(*,*) ' slave set:',slavset(1:iposs-1) 2493 write(*,*) ' master set:',mastset(1:iposm-1) 2494 call inputerror(inpc,ipoinpc,iline, 2495 & "*TIE%",ier) 2496 cycle loop 2497 endif 2498 call getnewline(inpc,textpart,istat,n,key,iline,ipol,inl, 2499 & ipoinp,inp,ipoinpc) 2500 elseif(textpart(1)(1:11).eq.'*TIMEPOINTS') then 2501 igen=.false. 2502 nam_=nam_+1 2503 do i=2,n 2504 if(textpart(i)(1:8).eq.'GENERATE') then 2505 igen=.true. 2506 exit 2507 endif 2508 enddo 2509 do 2510 call getnewline(inpc,textpart,istat,n,key,iline,ipol,inl, 2511 & ipoinp,inp,ipoinpc) 2512 if((istat.lt.0).or.(key.eq.1)) exit 2513 if(igen)then 2514 if(n.lt.3)then 2515 write(*,*)'*ERROR in allocation:' 2516 call inputerror(inpc,ipoinpc,iline, 2517 & "*TIME POINTS%",ier) 2518 exit 2519 else 2520 read(textpart(1)(1:20),'(f20.0)',iostat=istat) 2521 & tpmin 2522 if(istat.gt.0) then 2523 call inputerror(inpc,ipoinpc,iline, 2524 & "*TIME POINTS%",ier) 2525 exit 2526 endif 2527 read(textpart(2)(1:20),'(f20.0)',iostat=istat) 2528 & tpmax 2529 if(istat.gt.0) then 2530 call inputerror(inpc,ipoinpc,iline, 2531 & "*TIME POINTS%",ier) 2532 exit 2533 endif 2534 read(textpart(3)(1:20),'(f20.0)',iostat=istat) 2535 & tpinc 2536 if(istat.gt.0) then 2537 call inputerror(inpc,ipoinpc,iline, 2538 & "*TIME POINTS%",ier) 2539 exit 2540 endif 2541! 2542 if((tpinc.le.0).or.(tpmin.ge.tpmax)) then 2543 write(*,*) '*ERROR in allocation:' 2544 call inputerror(inpc,ipoinpc,iline, 2545 & "*TIME POINTS%",ier) 2546 exit 2547 else 2548 namtot_=namtot_+2+INT((tpmax-tpmin)/tpinc) 2549 endif 2550 2551 endif 2552 else 2553 namtot_=namtot_+8 2554 endif 2555 enddo 2556 elseif(textpart(1)(1:10).eq.'*TRANSFORM') then 2557 ntrans_=ntrans_+1 2558 do 2559 call getnewline(inpc,textpart,istat,n,key,iline,ipol,inl, 2560 & ipoinp,inp,ipoinpc) 2561 if((istat.lt.0).or.(key.eq.1)) exit 2562 enddo 2563 elseif(textpart(1)(1:11).eq.'*TRANSFORMF') then 2564 ntrans_=ntrans_+1 2565 surface(1:1)=' ' 2566 do i=2,n 2567 if(textpart(i)(1:8).eq.'SURFACE=') then 2568 surface=textpart(i)(9:88) 2569 ipos=index(surface,' ') 2570 surface(ipos:ipos)='T' 2571 exit 2572 endif 2573 enddo 2574 if(surface(1:1).ne.' ') then 2575c do i=1,nset_ 2576c if(set(i).eq.surface) then 2577 call cident81(set,surface,nset_,i) 2578 if(i.gt.0) then 2579 if(set(i).eq.surface) then 2580 nload_=nload_+meminset(i) 2581c exit 2582 endif 2583 endif 2584 endif 2585 do 2586 call getnewline(inpc,textpart,istat,n,key,iline,ipol,inl, 2587 & ipoinp,inp,ipoinpc) 2588 if((istat.lt.0).or.(key.eq.1)) exit 2589 enddo 2590 elseif(textpart(1)(1:12).eq.'*USERELEMENT') then 2591 call userelements(textpart,n,iuel,nuel,inpc,ipoinpc,iline, 2592 & ier,ipoinp,inp,inl,ipol) 2593 call getnewline(inpc,textpart,istat,n,key,iline,ipol,inl, 2594 & ipoinp,inp,ipoinpc) 2595 elseif(textpart(1)(1:13).eq.'*USERMATERIAL') then 2596 ntmatl=0 2597 do i=2,n 2598 if(textpart(i)(1:10).eq.'CONSTANTS=') then 2599 read(textpart(i)(11:20),'(i10)',iostat=istat) 2600 & nconstants 2601 if(istat.gt.0) then 2602 call inputerror(inpc,ipoinpc,iline, 2603 & "*USER MATERIAL%",ier) 2604 exit 2605 endif 2606 ncmat_=max(nconstants,ncmat_) 2607 exit 2608 endif 2609 enddo 2610 do 2611 call getnewline(inpc,textpart,istat,n,key,iline,ipol,inl, 2612 & ipoinp,inp,ipoinpc) 2613 if((istat.lt.0).or.(key.eq.1)) exit 2614 ntmatl=ntmatl+1 2615 ntmat_=max(ntmatl,ntmat_) 2616 do i=2,nconstants/8+1 2617 call getnewline(inpc,textpart,istat,n,key,iline,ipol, 2618 & inl,ipoinp,inp,ipoinpc) 2619 enddo 2620 enddo 2621 elseif(textpart(1)(1:12).eq.'*USERSECTION') then 2622 nconstants=0 2623 do i=2,n 2624 if(textpart(i)(1:10).eq.'CONSTANTS=') then 2625 read(textpart(i)(11:20),'(i10)',iostat=istat) 2626 & nconstants 2627 if(istat.gt.0) then 2628 call inputerror(inpc,ipoinpc,iline, 2629 & "*USER SECTION%",ier) 2630 exit 2631 endif 2632 nprop_=nprop_+nconstants 2633 exit 2634 endif 2635 enddo 2636 do 2637 call getnewline(inpc,textpart,istat,n,key,iline,ipol,inl, 2638 & ipoinp,inp,ipoinpc) 2639 if((istat.lt.0).or.(key.eq.1)) exit 2640 enddo 2641 else 2642! 2643 call getnewline(inpc,textpart,istat,n,key,iline,ipol,inl, 2644 & ipoinp,inp,ipoinpc) 2645 endif 2646 enddo loop 2647! 2648 do i=1,nset_ 2649 nalset_=nalset_+rmeminset(i) 2650 maxrmeminset=max(maxrmeminset,rmeminset(i)) 2651 enddo 2652! 2653! extra space needed for rearrangement in elements.f and 2654! noelsets.f 2655! 2656 nalset_=nalset_+maxrmeminset 2657! 2658 nmpc_=nmpc_+1 2659 memmpc_=memmpc_+1 2660! 2661 if(irstrt(1).eq.0) then 2662 ne1d=neb32 2663 ne2d=necper+necpsr+necaxr+nesr 2664 endif 2665! 2666! introducing a fake tie for axisymmetric elements 2667! (needed for cavity radiation) 2668! 2669 if(necaxr.gt.0) ntie_=max(1,ntie_) 2670! 2671! providing space for the expansion of shell and beam elements 2672! to genuine volume elements (no distinction is made between 2673! linear and quadratic elements. The worst case (quadratic) 2674! is taken 2675! 2676 nk_=nk_+3*8*ne2d+8*3*ne1d 2677 if(ne1d.gt.0) then 2678 nboun_=nboun_*9 2679 nforc_=nforc_*9 2680 elseif(ne2d.gt.0) then 2681 nboun_=nboun_*4 2682 nforc_=nforc_*4 2683 endif 2684! 2685! providing for rigid nodes (knots) 2686! 2687! number of knots: 8*ne2d+3*ne1d 2688! number of expanded nodes: 3*8*ne2d+8*3*ne1d 2689! 2690! number of extra nodes (1 rotational node and 2691! 1 expansion node per knot 2692! and one inhomogeneous term node per expanded node) 2693! 2694 nk_=nk_+(2+3)*8*ne2d+(2+8)*3*ne1d 2695! 2696! number of equations (3 per expanded node) 2697! 2698 nmpc_=nmpc_+3*(3*8*ne2d+8*3*ne1d) 2699! 2700! number of terms: 9 per equation 2701! 2702 memmpc_=memmpc_+9*3*(3*8*ne2d+8*3*ne1d) 2703! 2704! number of SPC's: 1 per DOF per expanded node 2705! 2706 nboun_=nboun_+3*(3*8*ne2d+8*3*ne1d) 2707! 2708! temperature DOF in knots 2709! 2710 nmpc_=nmpc_+(3*8*ne2d+8*3*ne1d) 2711 memmpc_=memmpc_+2*(3*8*ne2d+8*3*ne1d) 2712! 2713! extra MPCs to avoid undefinid rotation of rigid body nodes 2714! lying on a line 2715! 2716 nmpc_=nmpc_+3*8*ne2d+8*3*ne1d 2717 memmpc_=memmpc_+3*(3*8*ne2d+8*3*ne1d) 2718! 2719! extra nodes for the radiation boundary conditions 2720! 2721 nk_=nk_+nradiate 2722! 2723! each layer in each shell has a local orientation 2724! 2725 norien_=norien_+nesr*mi(3) 2726! 2727 write(*,*) 2728 write(*,*) ' The numbers below are estimated upper bounds' 2729 write(*,*) 2730 write(*,*) ' number of:' 2731 write(*,*) 2732 write(*,*) ' nodes: ',nk_ 2733 write(*,*) ' elements: ',ne_ 2734 write(*,*) ' one-dimensional elements: ',ne1d 2735 write(*,*) ' two-dimensional elements: ',ne2d 2736 write(*,*) ' integration points per element: ',mi(1) 2737 write(*,*) ' degrees of freedom per node: ',mi(2) 2738 write(*,*) ' layers per element: ',mi(3) 2739 write(*,*) 2740 write(*,*) ' distributed facial loads: ',nload_ 2741 write(*,*) ' distributed volumetric loads: ',nbody_ 2742 write(*,*) ' concentrated loads: ',nforc_ 2743 write(*,*) ' single point constraints: ',nboun_ 2744 write(*,*) ' multiple point constraints: ',nmpc_ 2745 write(*,*) ' terms in all multiple point constraints: ',memmpc_ 2746 write(*,*) ' tie constraints: ',ntie_ 2747 write(*,*) ' dependent nodes tied by cyclic constraints: ',ncs_ 2748 write(*,*) ' dependent nodes in pre-tension constraints: ',npt_ 2749 write(*,*) 2750 write(*,*) ' sets: ',nset_ 2751 write(*,*) ' terms in all sets: ',nalset_ 2752 write(*,*) 2753 write(*,*) ' materials: ',nmat_ 2754 write(*,*) ' constants per material and temperature: ',ncmat_ 2755 write(*,*) ' temperature points per material: ',ntmat_ 2756 write(*,*) ' plastic data points per material: ',npmat_ 2757 write(*,*) 2758 write(*,*) ' orientations: ',norien_ 2759 write(*,*) ' amplitudes: ',nam_ 2760 write(*,*) ' data points in all amplitudes: ',namtot_ 2761 write(*,*) ' print requests: ',nprint_ 2762 write(*,*) ' transformations: ',ntrans_ 2763 write(*,*) ' property cards: ',nprop_ 2764 write(*,*) 2765! 2766 if(ier.eq.1) then 2767 write(*,*) '*ERROR in allocation: at least one fatal' 2768 write(*,*) ' error message while reading the' 2769 write(*,*) ' input deck: CalculiX stops.' 2770 write(*,*) 2771 call exit(201) 2772 endif 2773! 2774 return 2775 end 2776