c   This program performs power calculations
c   for a proposed ASTM plastic lumber standard.
c   Steve Verrill, 3/25/99

      real xmoduli(5),cv(5),probf(5,5)

      integer icv(5),moduli(5)

      real*8 p,q,tval,pf,dfn,dfd,xncp,bound

      integer istatus

      character*80 results

      character*80 argv,intwrite

      call getarg(6,results)


      call getarg(7,results)

      write(6,9900) results
 9900 format("The results from this run have been appended to",
     x       " the file ",a80,"
To see how to download", x " this file,", x ' read this description of anonymous ftp', x '.

','Alternatively, you should be able to do a', x ' "Save As" or "Print" in your browser.

') write(7,9700) 9700 format(////,1x,"**************************************",//,1x, x "Output from astmplast1.f",//,1x, x "**************************************",/) c c secant moduli c call getarg(1,argv) write(intwrite,7) argv 7 format(a80) read(intwrite,*) nummod call getarg(2,argv) write(intwrite,7) argv read(intwrite,*) (moduli(i),i = 1,nummod) do 40 i = 1,nummod xmoduli(i) = moduli(i) 40 continue c c cvs c call getarg(3,argv) write(intwrite,7) argv read(intwrite,*) numcv call getarg(4,argv) write(intwrite,7) argv read(intwrite,*) (icv(i),i = 1,numcv) do 55 i = 1,numcv cv(i) = abs(icv(i)/100.0) if ((cv(i) .eq. 0.0) .or. (cv(i) .ge. 1.0)) then write(6,54) 54 format("
A coefficient of variation value", x " did not lie between 0 and 100. Please correct", x " the form.
") stop endif 55 continue c c number of specimens c call getarg(5,argv) write(intwrite,7) argv read(intwrite,*) num idof = num - 1 xn = num write(6,61) num write(7,62) num 61 format('
These calculations assume that the secant moduli', x' are normally distributed, and that ',i3,' specimens', x' will be tested.
') 62 format(/,1x,'These calculations assume that the secant moduli', x' are normally distributed,',/,1x, x'and that ',i3,' specimens will be tested.',/) c c perform the power calculations c fac = sqrt(xn) dfn = 1.0d0 dfd = xn - 1.0d0 p = .75d0 q = 1.0d0 - p call cdft(2,p,q,tval,dfd,istatus,bound) if (istatus .ne. 0) then write(6,71) istatus write(7,72) istatus 71 format('
Error 71: istatus = ',i2,'.
', x 'Please contact Steve Verrill at', x '', x 'sverrill@fs.fed.us.
') 72 format(/,1x,'Error 71: istatus = ',i2,'.',/, x 1x,'Please contact Steve Verrill at', x ' steve@swst.org.',/) stop endif do 110 i = 1,numcv do 100 j = 1,nummod xncp = fac*(xmoduli(j) - 50000)/ x (xmoduli(j)*cv(i)) c write(6,90) i,j,fac,xmoduli(j),cv(i),xncp c 90 format("
The i,j,fac,xmoduli(j),cv(i),xncp", c x " values are

",i1,2x,i1,2x,4(e12.4,2x), c x "
") ier = 0 call mdtn(tval,idof,xncp,pf,ier) if (ier .ne. 0) then write(6,81) istatus write(7,82) istatus 81 format('
Error 81: istatus = ',i2,'.
', x 'Please contact Steve Verrill at', x '', x 'sverrill@fs.fed.us.
') 82 format(/,1x,'Error 81: istatus = ',i2,'.',/, x 1x,'Please contact Steve Verrill at', x ' steve@swst.org.',/) stop endif probf(i,j) = pf 100 continue 110 continue c write(6,115) c 115 format("
Got to 115.
") write(6,1010) write(7,1012) (moduli(i), i = 1,nummod) 1010 format('

Moduli are on top. Coefficients', x' of variation are along the left side. Probabilities', x' of failing the test are in the table.

') 1012 format(///,1x,'Moduli are on top. Coefficients of', x' variation',/,1x,'are along the left side. ', x'Probabilities of failing the test',/,1x,'are in the table.', x//, x5x,5(2x,i6),/) do 1016 i = 1,numcv write(7,1014) icv(i),(probf(i,j),j = 1,nummod) 1014 format(2x,i3,5(3x,f5.3)) 1016 continue write(6,1020) 1020 format("") write(6,1030) nummod 1030 format("", x"") do 1034 i = 1,nummod write(6,1032) moduli(i) 1032 format("") 1034 continue write(6,1036) 1036 format("") do 1050 i = 1,numcv write(6,1042) icv(i),(probf(i,j), j = 1,nummod) 1042 format("","",("")) write(6,1044) 1044 format("") 1050 continue write(6,1052) 1052 format("
") write(6,6060) 6060 format('

For further information,', x' please contact Steve Verrill at ', x'', x'sverrill@fs.fed.us', x' or 608-231-9375.
Last modified on 3/25/99.') stop end