c ***************************************************************** c * * c * p o s t p r o c e s s o r p r o g r a m * c * * c ***************************************************************** implicit real*8 (a-h,o-z) common z(1000000) common /iotp/ nr,nw,lcp,lmd,lnq c----------------------------------------------------------------------- c This routine drives the postprocessing of the finite element analysis c----------------------------------------------------------------------- c c-----Set file numbers c ---------------- call setfil c c-----Read control parameters c ----------------------- call readcp c c-----Read output requests c -------------------- call reqinp (nrqtot) c c-----Initialize postprocessing c ------------------------- jktim = jdbadr ('ktim',nw) call initia (nbstp,z(jktim)) c c-----Loop over the time steps c ======================== nbarch = 0 nbrep = 0 ifirst = 1 do 100 istep = 1 , nbstp call timstp (ifirst,nbrep,nbarch,itnum,nrqtot) if (itnum .gt. 0) go to 150 100 continue c 150 continue c stop end c======================================================================= subroutine setfil c======================================================================= implicit real*8 (a-h,o-z) common /iotp/ nr,nw,lcp,lmd,lnq c----------------------------------------------------------------------- c This routine opens the external files needed for the program c----------------------------------------------------------------------- c c-----set files c --------- nr = 1 nw = 2 lcp = 7 lnq = 8 lmd = 9 open (nr ,file='pst.ips',status='old') open (nw ,file='pst.ops',status='new') open (lcp,file='prf.lcp',status='old',form='unformatted') open (lmd,file='fem.lmd',status='old',form='unformatted') open (lnq,file='pst.lnq',status='new') c return end c======================================================================= subroutine header (nr,nw,secnam) c======================================================================= implicit real*8 (a-h,o-z) character*6 atsign,secnam character*12 head,headr data atsign/'@@@@@@'/ c----------------------------------------------------------------------- c This routine read the standard atsign header for each new section c----------------------------------------------------------------------- c c-----Read header c ----------- read (nr,2000,err=900) head c c-----Check section header c -------------------- headr = atsign // secnam if (head .ne. headr) then write (nw,1000) headr,head stop endif c return c 900 write (nw,1001) secnam stop c 1000 format (///5x,'$$$ ERROR: INCONSISTENT DATA ENTRY DETECTED:'/ & 5x,'$$$ ',a12,' HEADER IS EXPECTED; ',a12,' WAS READ.'/ & 5x,'$$$ CHECK INPUT DATA SEQUENCE') 1001 format (///5x, & '$$$ AN READ ERROR OCCURED WHILE HEADING FOR SECTION ',a6) 2000 format (a12) end c==================================================================== subroutine rerror (secnam) c==================================================================== implicit real*8 (a-h,o-z) character*6 secnam common /iotp/ nr,nw,lcp,lmd,lnq c-------------------------------------------------------------------- c This routine processes read error messages. c-------------------------------------------------------------------- c c-----print out error message c ======================= write (nw,1000) secnam stop c 1000 format (///5x, & '$$$ AN INPUT ERROR OCCURED WHILE READING SECTION ',a6,' $$$') end c======================================================================= subroutine chkvar (varn,ivar,minvar,maxvar,iflag) c======================================================================= implicit real*8 (a-h,o-z) character*6 varn common /iotp/ nr,nw,lcp,lmd,lnq c----------------------------------------------------------------------- c This routine checks the bounds for integer variables. c----------------------------------------------------------------------- c c-----Check lower bound c ----------------- if (ivar .lt. minvar) then ivar = minvar write (nw,1000) varn,minvar if (iflag .gt. 0) stop endif c c-----Check upper bound c ----------------- if (ivar .gt. maxvar) then ivar = maxvar write (nw,1001) varn,maxvar if (iflag .gt. 0) stop endif c return c 1000 format (///10x, & '$$$ VARIABLE _',a6,'_ IS BELOW ALLOWABLE LIMIT $$$'/10x, & '$$$ IT WAS SET TO LOWER LIMIT =',i10,' $$$'///) 1001 format (///10x, & '$$$ VARIABLE _',a6,'_ IS ABOVE ALLOWABLE LIMIT $$$'/10x, & '$$$ IT WAS SET TO UPPER LIMIT =',i10,' $$$'///) end c======================================================================= subroutine readcp c======================================================================= implicit real*8 (a-h,o-z) dimension title(12) common /iotp/ nr,nw,lcp,lmd,lnq common /strc/ method,nnod,ntriad,nbody,nbvel,noddof, & ndof,ndedld,scalef common /aero/ nrotor,nlftln,ninfst,ntastn,idstal common /junk/ zero,half,one,two,four,pi,rmchep,small c----------------------------------------------------------------------- c This routine reads in the problem control parameters c----------------------------------------------------------------------- c c-----Read title c ---------- call header (nr,nw,' PST-1') read (nr,2000,err=900) title write (nw,1000 ) title c c-----Read control parameters c ----------------------- read (lcp) method,nnod,ntriad,nbody,nbvel,noddof, & ndof,ndedld,scalef read (lcp) nrotor,nlftln,ninfst,ntastn,idstal read (lcp) zero,half,one,two,four,pi,rmchep,small c c-----Retrieve data base c ------------------ call dbsret (lcp) c c-----Start reading archival file c --------------------------- read (lmd) ii c return c 900 call rerror ('PST-1 ') c 1000 format ( & 15x,'**************************************************'/ & 15x,'* CENTER OF EXCELLENCE IN ROTARY WING TECHNOLOGY *'/ & 15x,'* RENSSELAER POLYTECHNIC INSTITUTE *'/ & 15x,'* *'/ & 15x,'* - - - D Y M O R E - - - *'/ & 15x,'* *'/ & 15x,'* POSTPROCESSING OF THE FINITE ELEMENT ANALYSIS *'/ & 15x,'**************************************************'// & 5x,'{PST-1}: ',12a6) 2000 format (12a6) end c======================================================================= subroutine reqinp (nrqtot) c======================================================================= implicit real*8 (a-h,o-z) common z(1000000) common /iotp/ nr,nw,lcp,lmd,lnq common /strc/ method,nnod,ntriad,nbody,nbvel,noddof, & ndof,ndedld,scalef common /aero/ nrotor,nlftln,ninfst,ntastn,idstal c----------------------------------------------------------------------- c This routine drives the reading of the output requests c----------------------------------------------------------------------- c c-----Process storage requirements c ---------------------------- jnbel = jdbadr ('nbel',nw) jnbcl = jdbadr ('nbcl',nw) call dbsreq ('nrqc', 8,jnrqc,0,nw) call dbsreq ('nrql', 12,jnrql,0,nw) call dbsreq ('nrqd', 12,jnrqd,0,nw) call dbsreq ('kocp',nbody*4,jkocp,0,nw) call dbsreq ('iprs', 6,jiprs,0,nw) c c-----Read output control parameters c ============================== call outctp (nnod,ntriad,nbody,iprtbl,z(jkocp)) c c-----Read number of output requests c ------------------------------ call reqnbr (nnod,nbody,nrqtot,z(jnbel),z(jnbcl),z(jnrqc), & z(jnrql),z(jnrqd)) c c-----Read configuration output requests c ---------------------------------- ireqst = 0 call reqcnf (nnod,nlftln,ireqst,z(jnrqc)) c c-----Read output requests for structural elements c -------------------------------------------- call reqstr (ireqst,z(jnbel),z(jnrql)) c c-----Read output requests for constraint elements c -------------------------------------------- call reqcns (ireqst,z(jnbcl),z(jnrqd)) c return end c======================================================================= subroutine outctp (nnod,ntriad,nbody,iprtbl,kocp) c======================================================================= implicit real*8 (a-h,o-z) dimension kocp(1) common /iotp/ nr,nw,lcp,lmd,lnq c----------------------------------------------------------------------- c This routine read the number of output requests c----------------------------------------------------------------------- call header (nr,nw,' PST-2') write (nw,1000) c c-----Read output control parameters for each body c -------------------------------------------- do 10 ibody = 1 , nbody read (nr,*,err=900) k,irtrep,itrref,nodref write (nw, 1001) ibody,irtrep,itrref,nodref c c-----Check validity of data c ---------------------- call chkvar ('irtrep',irtrep,0, 1,0) call chkvar ('itrref',itrref,0,ntriad,0) call chkvar ('nodref',nodref,0, nnod,0) c c-----Save output control parameters c ------------------------------ ii = (ibody - 1) * 4 kocp(ii+1) = irtrep kocp(ii+2) = itrref kocp(ii+3) = nodref c 10 continue c return c 900 call rerror ('PST-2 ') c 1000 format (///20x,'***************************************'/ & 20x,'* OUTPUT CONTROL PARAMETERS {PST-2} *'/ & 20x,'***************************************'/// & ' BODY ROTATION REPRESENTATION REFERENCE REFERENCE'/ & ' NUMBER FLAG TRIAD NODE') 1001 format (i7,i27,2i12) end c======================================================================= subroutine reqnbr (nnod,nbody,nrqtot,nbel,nbcl,nrqc,nrql,nrqd) c======================================================================= implicit real*8 (a-h,o-z) dimension nbel(1),nbcl(1),nrqc(1),nrql(1),nrqd(1) common /iotp/ nr,nw,lcp,lmd,lnq c----------------------------------------------------------------------- c This routine read the number of output requests c----------------------------------------------------------------------- c c-----Read number of output requests for configuration c ================================================ call header (nr,nw,' REQ-1') read (nr,*,err=900) (nrqc(i),i=1,4) write (nw,1000) (nrqc(i),i=1,4) c c-----Check validity of data c ---------------------- call chkvar ('nrqdis',nrqc(1),0, nnod*6,0) call chkvar ('nrqerg',nrqc(2),0, 9,0) c c-----Read output request number for structural elements c ================================================== read (nr,*,err=900) (nrql(i),i=1,4) write (nw,1001) (nrql(i),i=1,4) c c-----Check validity of data c ---------------------- call chkvar ('nrqbdl',nrql( 1),0,nbel( 1)*64,0) call chkvar ('nrqrgb',nrql( 3),0,nbel( 3)* 6,0) call chkvar ('nrqflx',nrql( 4),0,nbel( 4)* 6,0) c c-----Read output request number for constraint elements c ================================================== read (nr,*,err=900) (nrqd(i),i=1,8) write (nw,1002) (nrqd(i),i=1,7) c c-----Check validity of data c ---------------------- call chkvar ('nrqprd',nrqd( 1),0,nbcl( 1)* 2,0) call chkvar ('nrqlcn',nrqd( 2),0,nbcl( 2)* 6,0) call chkvar ('nrqrvj',nrqd( 3),0,nbcl( 3)* 4,0) call chkvar ('nrqunj',nrqd( 4),0,nbcl( 4)* 4,0) call chkvar ('nrqshj',nrqd( 5),0,nbcl( 5)* 4,0) call chkvar ('nrqprj',nrqd( 6),0,nbcl( 6)* 4,0) call chkvar ('nrqslj',nrqd( 7),0,nbcl( 7)* 4,0) call chkvar ('nrqlnk',nrqd( 8),0,nbcl( 8)* 4,0) c c-----Total number of requests c ======================== nrqtot = 0 do 10 i = 1 , 4 nrqtot = nrqtot + nrqc(i) 10 continue do 15 i = 1 , 4 nrqtot = nrqtot + nrql(i) 15 continue do 20 i = 1 , 8 nrqtot = nrqtot + nrqd(i) 20 continue write (nw,1003) nrqtot c c-----Process storage requirements c ---------------------------- call dbsreq ('reqt',nrqtot,jreqt,1,nw) c return c 900 call rerror ('REQ-1 ') c 1000 format (///20x,'***********************'/20x, & '* OUTPUT REQUESTS *'/20x,'***********************'/// & 5x,'--- NUMBER OF CONFIGURATION OUTPUT REQUESTS: {REQ-1}'// & 5x,'NODAL DISPLACEMENTS = ',i5/ & 5x,'ENERGIES = ',i5/ & 5x,'AIRLOADS = ',i5/ & 5x,'GENERALIZED INFLOW = ',i5) 1001 format (///5x, & '--- OUTPUT REQUEST NUMBERS FOR STRUCTURAL ELEMENTS: {REQ-1}'// & 5x,'BEAM ELEMENTS =',i5/ & 5x,'ADVANCED BEAM ELEMENTS =',i5/ & 5x,'RIGID BODIES =',i5/ & 5x,'FLEXIBLE JOINT ELEMENTS =',i5) 1002 format (///5x, & '--- OUTPUT REQUEST NUMBERS FOR CONSTRAINT ELEMENTS: {REQ-1}'// & 5x,'PRESCRIBED DISPLACEMENTS =',i5/ & 5x,'REVOLUTE JOINTS =',i5/ & 5x,'UNIVERSAL JOINTS =',i5/ & 5x,'SPHERICAL JOINTS =',i5/ & 5x,'PRISMATIC JOINTS =',i5/ & 5x,'SLIDING JOINTS =',i5/ & 5x,'RIGID LINK ELEMENT =',i5) 1003 format (//5x,'TOTAL NUMBER OF OUTPUT REQUESTS = ',i5) end c======================================================================= subroutine reqcnf (nnod,nlftln,ireqst,nrqc) c======================================================================= implicit real*8 (a-h,o-z) dimension nrqc(1) common z(1000000) common /iotp/ nr,nw,lcp,lmd,lnq c----------------------------------------------------------------------- c This routine reads the output requests for system configuration c----------------------------------------------------------------------- c c-----Read displacement output requests c --------------------------------- nrq = nrqc(1) call dbsreq ('iqds',nrq*5,jiqds,0,nw) if (nrq .gt. 0) then call reqdis (nnod,nrq,ireqst,z(jiqds)) endif c c-----Read energy output requests c --------------------------- nrq = nrqc(2) call dbsreq ('iqrg',nrq*2,jiqrg,0,nw) if (nrq .gt. 0) then call reqerg (nrq,ireqst,z(jiqrg)) endif c c-----Read airstation output requests c ------------------------------- nrq = nrqc(3) call dbsreq ('iqas',nrq*3,jiqas,0,nw) if (nrq .gt. 0) then jklfn = jdbadr ('klfn',nw) call reqast (nlftln,nrq,ireqst,z(jklfn),z(jiqas)) endif c return end c======================================================================= subroutine reqdis (nnod,nrqdis,ireqst,iqds) c======================================================================= implicit real*8 (a-h,o-z) dimension iqds(1) common /iotp/ nr,nw,lcp,lmd,lnq c----------------------------------------------------------------------- c This routine reads in displacement output requests c----------------------------------------------------------------------- call header (nr,nw,' REQ-2') write (nw,1000) c c-----Read nodal displacement output requests c --------------------------------------- do 10 irqdis = 1 , nrqdis ireqst = ireqst + 1 read (nr,*,err=900) k,inod,idof,irtrep,iaxchg write (nw,1001) ireqst,inod,idof,irtrep,iaxchg c c-----Check validity of data c ---------------------- call chkvar (' inod', inod,1,nnod,0) if (irtrep .eq. 0) then call chkvar (' idof', idof,1, 6,0) else call chkvar (' idof', idof,1, 12,0) endif call chkvar ('irtrep',irtrep,0, 1,0) call chkvar ('iaxchg',iaxchg,0, 1,0) c ii = (irqdis - 1) * 5 iqds(ii+1) = ireqst iqds(ii+2) = inod iqds(ii+3) = idof iqds(ii+4) = irtrep iqds(ii+5) = iaxchg c 10 continue c return c 900 call rerror ('REQ-2 ') c 1000 format (///5x,'--- DISPLACEMENT OUTPUT REQUESTS {REQ-2}'// & ' REQUEST DISPLACEMENT DEGREE OF ROTATION AXIS'/ & ' NUMBER AT NODE NO FREEDOM REPRESENTATION SYSTEM'/) 1001 format (i8,i14,i11,i15,i7) end c======================================================================= subroutine reqerg (nrqerg,ireqst,iqrg) c======================================================================= implicit real*8 (a-h,o-z) dimension iqrg(1) common /iotp/ nr,nw,lcp,lmd,lnq c----------------------------------------------------------------------- c This routine reads in energy output requests c----------------------------------------------------------------------- call header (nr,nw,' REQ-3') write (nw,1000) c c-----Read energy output requests c --------------------------- do 10 irqerg = 1 , nrqerg ireqst = ireqst + 1 read (nr,*,err=900) k,idof c c-----Check validity of data c ---------------------- call chkvar (' idof', idof,1,9,0) write (nw,1001) ireqst,idof c c-----Write and save output request c ----------------------------- ii = (irqerg - 1) * 2 iqrg(ii+1) = ireqst iqrg(ii+2) = idof c 10 continue c return c 900 call rerror ('REQ-3 ') c 1000 format (///5x,'--- ENERGY OUTPUT REQUESTS {REQ-3}'// & ' REQUEST ENERGY'/ & ' NUMBER TYPE'/) 1001 format (i8,i8) end c======================================================================= subroutine reqast (nlftln,nrqast,ireqst,klfn,iqas) c======================================================================= implicit real*8 (a-h,o-z) dimension klfn(1),iqas(1) common /iotp/ nr,nw,lcp,lmd,lnq c----------------------------------------------------------------------- c this routine reads in the airloads output requests c----------------------------------------------------------------------- call header (nr,nw,' REQ-4') write (nw,1000) c c-----Read airloads output requests c ----------------------------- do 10 irqast = 1 , nrqast ireqst = ireqst + 1 read (nr,*,err=900) k,ilftln,iastn,idof c call chkvar ('ilftln',ilftln,1,nlftln,0) c c-----Extract aerodynamic reference line data c --------------------------------------- call lfnelp (ilftln,irotor,iairst,lairst,nodlfn,itrlfn,klfn) c c-----Check validity of data c ---------------------- ii = lairst - iairst + 1 call chkvar (' iastn',iastn,1,ii,0) call chkvar (' idof', idof,1,50,0) c c-----Write and save output request c ----------------------------- write (nw,1001) ireqst,ilftln,iastn,idof c ii = (irqast - 1) * 3 iqas(ii+1) = ireqst iqas(ii+2) = iairst + iastn - 1 iqas(ii+3) = idof c 10 continue c return c 900 call rerror (' REQ-4') c 1000 format (///5x,'--- AIRSTATION DATA OUTPUT REQUESTS {REQ-4}'// & ' REQUEST AERO REF. AIRSTATION DATA'/ & ' NUMBER LINE NUMBER NUMBER'/) 1001 format (i8,3i12) end c======================================================================= subroutine reqstr (ireqst,nbel,nrql) c======================================================================= implicit real*8 (a-h,o-z) dimension nbel(1),nrql(1) common z(1000000) common /iotp/ nr,nw,lcp,lmd,lnq c----------------------------------------------------------------------- c This routine reads the output requests for structural elements c----------------------------------------------------------------------- c c-----Internal quantities in blade elements c ------------------------------------- nbelm = nbel(1) nrqlm = nrql(1) call dbsreq ('iqbl',nrqlm*5,jiqbl,0,nw) if ( (nbelm .gt. 0) .and. (nrqlm .gt. 0) ) then call reqbld (nbelm,nrqlm,ireqst,z(jiqbl)) endif c c-----Internal quantities in rigid bodies c ----------------------------------- nbelm = nbel(3) nrqlm = nrql(3) call dbsreq ('iqrb',nrqlm*4,jiqrb,0,nw) if ( (nbelm .gt. 0) .and. (nrqlm .gt. 0) ) then call reqrgb (nbelm,nrqlm,ireqst,z(jiqrb)) endif c c-----Internal quantities in flexible joints c -------------------------------------- nbelm = nbel(4) nrqlm = nrql(4) call dbsreq ('iqfx',nrqlm*4,jiqfx,0,nw) if ( (nbelm .gt. 0) .and. (nrqlm .gt. 0) ) then call reqflx (nbelm,nrqlm,ireqst,z(jiqfx)) endif c return end c======================================================================= subroutine reqbld (nbld,nrqbld,ireqst,iqbl) c======================================================================= implicit real*8 (a-h,o-z) dimension iqbl(1) common /iotp/ nr,nw,lcp,lmd,lnq c----------------------------------------------------------------------- c This routine reads in the blade output requests c----------------------------------------------------------------------- call header (nr,nw,' RSE-1') write (nw,1000) c c-----Read internal quantities output requests c ---------------------------------------- do 10 irqbld = 1 , nrqbld ireqst = ireqst + 1 read (nr,*,err=900) k,noel,iosp,idof,iflg c c-----Check validity of data c ---------------------- call chkvar (' noel',noel,1,nbld,0) call chkvar (' iosp',iosp,1, 3,0) call chkvar (' idof',idof,1, 6,0) call chkvar (' iflg',iflg,1, 9,0) c write (nw,1001) ireqst,noel,iosp,idof,iflg c ii = (irqbld - 1) * 5 iqbl(ii+1) = ireqst iqbl(ii+2) = noel iqbl(ii+3) = iosp iqbl(ii+4) = idof iqbl(ii+5) = iflg c 10 continue c return c 900 call rerror ('RSE-1 ') c 1000 format (/// & 5x,'--- BLADE ELEMENT OUTPUT REQUESTS {RSE-1}'// & ' REQUEST BLADE ELEMENT SAMPLING DEGREE OF OPTION'/ & ' NUMBER NUMBER POINT NO FREEDOM FLAG'/) 1001 format (i8,i15,2i11,i8) end c======================================================================= subroutine reqrgb (nrgb,nrqrgb,ireqst,iqrb) c======================================================================= implicit real*8 (a-h,o-z) dimension iqrb(1) common /iotp/ nr,nw,lcp,lmd,lnq c----------------------------------------------------------------------- c This routine reads in the rigid body output requests c----------------------------------------------------------------------- call header (nr,nw,' RSE-3') write (nw,1000) c c-----Read internal quantities output requests c ---------------------------------------- do 10 irqrgb = 1 , nrqrgb ireqst = ireqst + 1 read (nr,*,err=900) k,noel,idof,iflg c c-----Check validity of data c ---------------------- call chkvar (' noel',noel,1,nrgb,0) call chkvar (' idof',idof,1, 6,0) call chkvar (' iflg',iflg,1, 3,0) c write (nw,1001) ireqst,noel,idof,iflg c ii = (irqrgb - 1) * 4 iqrb(ii+1) = ireqst iqrb(ii+2) = noel iqrb(ii+3) = idof iqrb(ii+4) = iflg c 10 continue c return c 900 call rerror ('RSE-3 ') c 1000 format (/// & 5x,'--- RIGID BODY OUTPUT REQUESTS {RSE-3}'// & ' REQUEST RIGID BODY DEGREE OF OPTION'/ & ' NUMBER NUMBER FREEDOM FLAG'/) 1001 format (i8,i13,i11,i8) end c======================================================================= subroutine reqflx (nflx,nrqflx,ireqst,iqfx) c======================================================================= implicit real*8 (a-h,o-z) dimension iqfx(1) common /iotp/ nr,nw,lcp,lmd,lnq c----------------------------------------------------------------------- c This routine reads in the flexible joint output requests c----------------------------------------------------------------------- call header (nr,nw,' RSE-4') write (nw,1000) c c-----Read internal quantities output requests c ---------------------------------------- do 10 irqflx = 1 , nrqflx ireqst = ireqst + 1 read (nr,*,err=900) k,noel,idof,iflg c c-----Check validity of data c ---------------------- call chkvar (' noel',noel,1,nflx,0) call chkvar (' idof',idof,1, 6,0) call chkvar (' iflg',iflg,1, 6,0) c write (nw,1001) ireqst,noel,idof,iflg c ii = (irqflx - 1) * 4 iqfx(ii+1) = ireqst iqfx(ii+2) = noel iqfx(ii+3) = idof iqfx(ii+4) = iflg c 10 continue c return c 900 call rerror ('RSE-4 ') c 1000 format (/// & 5x,'--- FLEXIBLE JOINT OUTPUT REQUESTS {RSE-4}'// & ' REQUEST FLEXIBLE JOINT DEGREE OF OPTION'/ & ' NUMBER NUMBER FREEDOM FLAG'/) 1001 format (i8,i16,i11,i8) end c======================================================================= subroutine reqcns (ireqst,nbcl,nrqd) c======================================================================= implicit real*8 (a-h,o-z) dimension nbcl(1),nrqd(1) common z(1000000) common /iotp/ nr,nw,lcp,lmd,lnq c----------------------------------------------------------------------- c This routine reads the output requests for internal quantities c----------------------------------------------------------------------- c c-----Internal quantities in prescribed displacements c ----------------------------------------------- nbclm = nbcl(1) nrqdm = nrqd(1) call dbsreq ('iqpd', nbclm*3,jiqpd,0,nw) if ( (nbclm .gt. 0) .and. (nrqdm .gt. 0) ) then call reqprd (nbclm,nrqdm,ireqst,z(jiqpd)) endif c c-----Internal quantities in linear constraints c ----------------------------------------- nbclm = nbcl(2) nrqdm = nrqd(2) call dbsreq ('iqlc', nbclm*3,jiqlc,0,nw) if ( (nbclm .gt. 0) .and. (nrqdm .gt. 0) ) then call reqlcn (nbclm,nrqdm,ireqst,z(jiqlc)) endif c c-----Internal quantities in revolute joints c -------------------------------------- nbclm = nbcl(3) nrqdm = nrqd(3) call dbsreq ('iqrv', nbclm*12,jiqrv,0,nw) if ( (nbclm .gt. 0) .and. (nrqdm .gt. 0) ) then call reqrvj (nbclm,nrqdm,ireqst,z(jiqrv)) endif c c-----Internal forces in prismatic joints c ----------------------------------- nbclm = nbcl(6) nrqdm = nrqd(6) call dbsreq ('iqpr',nrqdm*3,jiqpr,0,nw) if ( (nbclm .gt. 0) .and. (nrqdm .gt. 0) ) then call reqprj (nbclm,nrqdm,ireqst,z(jiqpr)) endif c c-----Internal quantities in link elements c ------------------------------------ nbclm = nbcl(8) nrqdm = nrqd(8) call dbsreq ('iqln',nrqdm*3,jiqln,0,nw) if ( (nbclm .gt. 0) .and. (nrqdm .gt. 0) ) then call reqlnk (nbclm,nrqdm,ireqst,z(jiqln)) endif c return end c======================================================================= subroutine reqprd (nprd,nrqprd,ireqst,iqpd) c======================================================================= implicit real*8 (a-h,o-z) dimension iqpd(1) common /iotp/ nr,nw,lcp,lmd,lnq c----------------------------------------------------------------------- c This routine reads in the prescribed displacement output requests c----------------------------------------------------------------------- call header (nr,nw,' RCE-1') write (nw,1000) c c-----Read internal torques output requests c ------------------------------------- do 10 irqprd = 1 , nrqprd ireqst = ireqst + 1 read (nr,*,err=900) k,noel,idof c c-----Check validity of data c ---------------------- call chkvar (' noel',noel,1,nprd,0) call chkvar (' idof',idof,1, 2,0) c write (nw,1001) ireqst,noel,idof c ii = (irqprd - 1) * 3 iqpd(ii+1) = ireqst iqpd(ii+2) = noel iqpd(ii+3) = idof c 10 continue c return c 900 call rerror ('RCE-1 ') c 1000 format (///5x, & '--- PRESCRIBED DISPLACEMENT OUTPUT REQUESTS {RCE-1}'/// & ' REQUEST PRESCRIBED DISPLACEMENT DEGREE OF'/ & ' NUMBER NUMBER FREEDOM'/) 1001 format (i8,i25,i11) end c======================================================================= subroutine reqlcn (nlcn,nrqlcn,ireqst,iqlc) c======================================================================= implicit real*8 (a-h,o-z) dimension iqlc(1) common /iotp/ nr,nw,lcp,lmd,lnq c----------------------------------------------------------------------- c This routine reads in the linear constraint output requests c----------------------------------------------------------------------- call header (nr,nw,' RCE-2') write (nw,1000) c c-----Read internal force output requests c ----------------------------------- do 10 irqlcn = 1 , nrqlcn ireqst = ireqst + 1 read (nr,*,err=900) k,noel,idof c c-----Check validity of data c ---------------------- call chkvar (' noel',noel,1,nlcn,0) call chkvar (' idof',idof,1, 6,0) c write (nw,1001) ireqst,noel,idof c ii = (irqlcn - 1) * 3 iqlc(ii+1) = ireqst iqlc(ii+2) = noel iqlc(ii+3) = idof c 10 continue c return c 900 call rerror ('RCE-2 ') c 1000 format (///5x,'--- LINEAR CONSTRAINT OUTPUT REQUESTS {RCE-2}'/// & ' REQUEST LINEAR CONSTRAINT'/ & ' NUMBER NUMBER'/) 1001 format (i8,i19,i11) end c======================================================================= subroutine reqrvj (nrvj,nrqrvj,ireqst,iqrv) c======================================================================= implicit real*8 (a-h,o-z) dimension iqrv(1) common /iotp/ nr,nw,lcp,lmd,lnq c----------------------------------------------------------------------- c This routine reads in the revolute joint torque output requests c----------------------------------------------------------------------- call header (nr,nw,' RCE-3') write (nw,1000) c c-----Read internal torques output requests c ------------------------------------- do 10 irqrvj = 1 , nrqrvj ireqst = ireqst + 1 read (nr,*,err=900) k,noel,idof c c-----Check validity of data c ---------------------- call chkvar (' noel',noel,1,nrvj,0) call chkvar (' idof',idof,1, 4,0) c write (nw,1001) ireqst,noel,idof c ii = (irqrvj - 1) * 3 iqrv(ii+1) = ireqst iqrv(ii+2) = noel iqrv(ii+3) = idof c 10 continue c return c 900 call rerror ('RCE-3 ') c 1000 format (///5x,'--- REVOLUTE JOINT OUTPUT REQUESTS {RCE-3}'/// & ' REQUEST REVOLUTE JOINT DEGREE OF'/ & ' NUMBER NUMBER FREEDOM'/) 1001 format (i8,i16,i11) end c======================================================================= subroutine reqprj (nprj,nrqprj,ireqst,iqp) c======================================================================= implicit real*8 (a-h,o-z) dimension iqp(1) common /iotp/ nr,nw,lcp,lmd,lnq c----------------------------------------------------------------------- c This routine reads in the prismatic joints output requests c----------------------------------------------------------------------- call header (nr,nw,' RCE-6') write (nw,1000) c c-----Read internal loads output requests c ----------------------------------- do 10 irqprj = 1 , nrqprj ireqst = ireqst + 1 read (nr,*,err=900) k,noel,idof c c-----Check validity of data c ---------------------- call chkvar (' noel',noel,1,nprj,0) call chkvar (' idof',idof,1, 7,0) c write (nw,1001) ireqst,noel,idof c ii = (irqprj - 1) * 3 iqp(ii+1) = ireqst iqp(ii+2) = noel iqp(ii+3) = idof c 10 continue c return c 900 call rerror ('RCE-6 ') c 1000 format (///5x,'--- PRISMATIC JOINT OUTPUT REQUESTS {RCE-6}'/// & ' REQUEST PRISMATC JOINT DEGREE OF'/ & ' NUMBER NUMBER FREEDOM'/) 1001 format (i8,i16,i11) end c======================================================================= subroutine reqlnk (nlnk,nrqlnk,ireqst,iqln) c======================================================================= implicit real*8 (a-h,o-z) dimension iqln(1) common /iotp/ nr,nw,lcp,lmd,lnq c----------------------------------------------------------------------- c This routine reads in the link load output requests c----------------------------------------------------------------------- call header (nr,nw,' RCE-7') write (nw,1000) c c-----Read internal forces output requests c ------------------------------------ do 10 irqlnk = 1 , nrqlnk ireqst = ireqst + 1 read (nr,*,err=900) noel,idof c c-----Check validity of data c ---------------------- call chkvar (' noel',noel,1,nlnk,0) call chkvar (' idof',idof,1, 3,0) c write (nw,1001) ireqst,noel,idof c ii = (irqlnk - 1) * 3 iqln(ii+1) = ireqst iqln(ii+2) = noel iqln(ii+3) = idof c 10 continue c return c 900 call rerror ('RCE-7 ') c 1000 format (///5x,'--- LINK LOAD OUTPUT REQUESTS {RCE-7}'// & ' REQUEST LINK ELEMENT DEGREE OF'/ & ' NUMBER NUMBER FREEDOM'/) 1001 format (i8,i15,i11) end c======================================================================= subroutine initia (nbstp,ktim) c======================================================================= implicit real*8 (a-h,o-z) dimension ktim(1) common z(1000000) c----------------------------------------------------------------------- c This routine initializes the analysis c----------------------------------------------------------------------- c c-----Number of time steps c -------------------- nbstp = ktim(1) c c-----Initialize structural configuration c ----------------------------------- call inistr c return end c======================================================================= subroutine inistr c======================================================================= implicit real*8 (a-h,o-z) common z(1000000) common /iotp/ nr,nw,lcp,lmd,lnq common /strc/ method,nnod,ntriad,nbody,nbvel,noddof, & ndof,ndedld,scalef c----------------------------------------------------------------------- c This routine initializes the structural configuration c----------------------------------------------------------------------- c c-----Process structural storage requirements c --------------------------------------- nbdis = nnod*6 call dbsreq ('uun ',nbdis, juun,1,nw) call dbsreq ('duu ',nbdis, jduu,1,nw) call dbsreq ('veli',nbvel,jveli,1,nw) call dbsreq ('velf',nbvel,jvelf,1,nw) call dbsreq ('enrg', 9,jenrg,1,nw) c c-----Clear structural configuration c ------------------------------ call clrstr (nnod,nbvel,noddof,z(juun),z(jduu),z(jveli),z(jvelf), & z(jenrg)) c return end c======================================================================= subroutine clrstr (nnod,nbvel,noddof,uun,duu,veli,velf,enrg) c======================================================================= implicit real*8 (a-h,o-z) dimension uun(1),duu(1),veli(1),velf(1),enrg(1) common /junk/ zero,half,one,two,four,pi,rmchep,small c----------------------------------------------------------------------- c This routine initializes the structural configuration c----------------------------------------------------------------------- c c-----Initialize displacement and increments c -------------------------------------- nbdis = nnod * 6 do 10 i = 1 , nbdis uun(i) = zero 10 continue c nbdis = nnod * noddof do 12 i = 1 , nbdis duu(i) = zero 12 continue c c-----Initialize velocities c --------------------- do 20 i = 1 , nbvel veli(i) = zero velf(i) = zero 20 continue c c-----Initialize energies c ------------------- do 30 i = 1 , 9 enrg(i) = zero 30 continue c return end c======================================================================= subroutine timstp (ifirst,nbrep,nbarch,itnum,nrqtot) c======================================================================= implicit real*8 (a-h,o-z) common z(1000000) common /iotp/ nr,nw,lcp,lmd,lnq c----------------------------------------------------------------------- c This routine drives the time stepping analysis c----------------------------------------------------------------------- jnrqc = jdbadr ('nrqc',nw) jreqt = jdbadr ('reqt',nw) jiprs = jdbadr ('iprs',nw) jduu = jdbadr ('duu ',nw) c c-----Read control parameters c ----------------------- call stpinp (istep,nbrep,nbeig,iwrite,z(jiprs)) c c-----Post-process this time step c =========================== if (iwrite .eq. 0) then call ppsstp (ifirst,istep,itnum,nrqtot,z(jiprs),z(jnrqc), & z(jreqt)) endif c c-----Post-process eigenvectors c ------------------------- if (nbeig .gt. 0) then call eigpro (nbeig,z(jduu)) endif c return end c======================================================================= subroutine stpinp (istep,nbrep,nbeig,iwrite,iprs) c======================================================================= implicit real*8 (a-h,o-z) dimension iprs(1) common /iotp/ nr,nw,lcp,lmd,lnq c----------------------------------------------------------------------- c This routine reads the control parameters for a time step c----------------------------------------------------------------------- c c-----Read new control parameters c =========================== if (nbrep .eq. 0) then call header (nr,nw,' PRN-1') read (nr,*,err=900) (iprs(i),i=1,6),nbrep c c-----Use old control parameters c ========================== else nbrep = nbrep - 1 endif c read (lmd,err=901) istep,nbeig,iwrite c return c 900 call rerror ('PRN-1 ') 901 write (nw,1000) stop c 1000 format (///5x,'$$$ END OF ARCHIVAL FILE WAS ENCOUNTERED'/ & 5x,'$$$ WHILE READING TIME STEP CONTROL PARAMETERS') end c======================================================================= subroutine ppsstp (ifirst,istep,itnum,nrqtot,iprs,nrqc,reqt) c======================================================================= implicit real*8 (a-h,o-z) dimension iprs(1),nrqc(1),reqt(1) common z(1000000) common /iotp/ nr,nw,lcp,lmd,lnq common /strc/ method,nnod,ntriad,nbody,nbvel,noddof, & ndof,ndedld,scalef common /aero/ nrotor,nlftln,ninfst,ntastn,idstal c----------------------------------------------------------------------- c This routine drives the printing for one time step. c----------------------------------------------------------------------- juun = jdbadr ('uun ',nw) jveli = jdbadr ('veli',nw) jvelf = jdbadr ('velf',nw) jenrg = jdbadr ('enrg',nw) jzast = jdbadr ('zast',nw) c c-----Read structural configuration c ============================= call stpstr (nnod,nbvel,itnum,timei,deltat,z(juun),z(jveli), & z(jvelf),z(jenrg)) timef = timei + deltat c c-----Set initial energy c ------------------ if (ifirst .eq. 1) then call enrgin (ifirst,z(jenrg)) endif c c-----Read airstation data c ==================== if (ntastn .gt. 0) then call stpast (ntastn,z(jzast)) endif c c-----Print out time step heading c =========================== if ( (iprs(1) .gt. 0) .or. (iprs(2) .gt. 0) .or. & (iprs(3) .gt. 0) .or. (iprs(4) .gt. 0) ) then write (nw,1000) istep,timef if (method .eq. 0) then write (nw,1001) iabs(itnum) else write (nw,1002) iabs(itnum) endif endif c c-----Print out configuration c ======================= ieig = 0 call prncnf (nnod,nbody,noddof,ieig,scalef,deltat) c c-----Save output request table c ========================= if (nrqtot .gt. 0) then write (lnq,1003) timef,(reqt(i),i=1,nrqtot) endif c return c 1000 format(// & 23x,'*****************************************************'/ & 23x,'* TIME STEP NO: ',i5,' (AT TIME = ',1pe10.3,') *'/ & 23x,'*****************************************************') 1001 format (//5x,'STATIC EQUILIBRIUM CONFIGURATION ', & '(NB OF ITERATIONS = ',i3,')'/5x, & '================================') 1002 format (//5x,'DYNAMIC EQUILIBRIUM CONFIGURATION ', & '(NB OF ITERATIONS = ',i3,')'/5x, & '=================================') c1003 format (73(1x,1pe17.10)) 1003 format (73(1x,1pe23.16)) end c======================================================================= subroutine stpstr (nnod,nbvel,itnum,timei,deltat,uun,veli,velf, & enrg) c======================================================================= implicit real*8 (a-h,o-z) dimension uun(1),veli(1),velf(1),enrg(1) common /iotp/ nr,nw,lcp,lmd,lnq c----------------------------------------------------------------------- c This routine reads system structural configuration c----------------------------------------------------------------------- nbdis = nnod*6 c c-----Update velocities c ----------------- do 10 i = 1 , nbvel veli(i) = velf(i) 10 continue c c-----Readin configuration at this time c ================================= read (lmd,err=900) itnum,timei,deltat read (lmd,err=900) ( uun(i),i=1,nbdis) read (lmd,err=900) (velf(i),i=1,nbvel) read (lmd,err=900) (enrg(i),i=1, 8) c return c 900 write (nw,1000) stop c 1000 format (///5x,'$$$ END OF ARCHIVAL FILE WAS ENCOUNTERED'/ & 5x,'$$$ WHILE READING STRUCTURAL CONFIGURATION') end c======================================================================= subroutine stpast (ntastn,zast) c======================================================================= implicit real*8 (a-h,o-z) dimension zast(1) common /iotp/ nr,nw,lcp,lmd,lnq c----------------------------------------------------------------------- c This routine reads airstation data c----------------------------------------------------------------------- c c-----Readin airstation data c ---------------------- ii = ntastn*50 read (lmd,err=900) (zast(i),i=1,ii) c return c 900 write (nw,1000) stop c 1000 format (///5x,'$$$ END OF ARCHIVAL FILE WAS ENCOUNTERED'/ & 5x,'$$$ WHILE READING AIRSTATION DATA') end c======================================================================= subroutine enrgin (ifirst,enrg) c======================================================================= implicit real*8 (a-h,o-z) dimension enrg(1) c----------------------------------------------------------------------- c This routine initializes the total energy c----------------------------------------------------------------------- ifirst = 0 c c-----Initialize total energy c ----------------------- enrg(9) = enrg(1) + enrg(2) write (2,3000) (enrg(i),i=1,9) 3000 format ('enrg'/6(2x,1pe12.5)) c return end c======================================================================= subroutine eigpro (nbeig,eiv) c======================================================================= implicit real*8 (a-h,o-z) dimension eiv(1) common /iotp/ nr,nw,lcp,lmd,lnq common /strc/ method,nnod,ntriad,nbody,nbvel,noddof, & ndof,ndedld,scalef c----------------------------------------------------------------------- c This routine print out the eigensolution c----------------------------------------------------------------------- write (nw,1000) nbdis = nnod*6 c c-----Loop over eigenvectors c ====================== ieig = 1 500 if (ieig .le. nbeig) then c c-----Read eigenmode c -------------- read (lmd,err=900) are,aim read (lmd,err=900) (eiv(i),i=1,nbdis) write (nw,1001) ieig,are,aim c c-----Print real eigenmode c -------------------- deltat = zero call prncnf (nnod,nbody,noddof,ieig,scalef,deltat) c c-----Print imaginary eigenmode c ------------------------- if (aim .ne. zero) then write (nw,1002) ieig,are,aim read (lmd,err=900) (eiv(i),i=1,nbdis) call prncnf (nnod,nbody,noddof,ieig,scalef,deltat) ieig = ieig + 1 endif ieig = ieig + 1 c go to 500 endif c return c 900 write (nw,1003) stop c 1000 format (//5x, & 'EIGEN-ANALYSIS ABOUT THIS STATIC EQUILIBRIUM CONFIGURATION'/5x, & '==========================================================') 1001 format (//5x,'EIGENVALUE NUMBER:',i5,5x,'(',1pe12.5,',+-i ', & e12.5,') RAD/SEC; REAL PART OF EIGENMODE.'/5x, & '-----------------------') 1002 format (//5x,'EIGENVALUE NUMBER:',i5,5x,'('1pe12.5,',+-i ', & e12.5,') RAD/SEC; IMAGINARY PART OF EIGENMODE.'/5x, & '-----------------------') 1003 format (///5x,'$$$ END OF ARCHIVAL FILE WAS ENCOUNTERED'/ & 5x,'$$$ WHILE READING EIGENMODES') end