LVL PSL N
E.1 Simulation programs
The following sections introduce the four main simulation routines. The programs are fairly flexible in that the layup of LVL as well as the species mix used in PSL can be specified arbitrarily. Species-specific information is stored in the Species module (see section E.2) which can be expanded to include further species. Certain calculations and simulations are executed by calling the functions described in section E.3.
Important notice: the simulation routines as well as some of the functions use external functions from the statistical IMSL libraries. The name of these libraries must be specified in the project settings.
The LVLBend program
Figure E.1 shows the flowchart for the program LVLBend, which simulates the bending MOE of LVL. This flowchart, with the code included below clarifies the operation of the program almost fully. The array S contains the species layup from the compression to the tension side of the beam, coded by integer numbers. The species code is 1 for quaking aspen, 2 for red oak and 3 for yellow-poplar.
PROGRAM LVLBend USE Species USE Ifaces
CHARACTER (LEN=50) :: FReg,FEflat,FEedge,FThick,Fdens DOUBLE PRECISION X,W,T,Ifl,Ied,Efl,Eed,MP,A,NN
DOUBLE PRECISION, DIMENSION(3) :: Regpar DOUBLE PRECISION, DIMENSION(5) :: Stat
START
Layup specification
DO simulations
Generate for each layer:
Original thickness Layer thickness Layer weight
EC and ET (including the densification effect) INPUT random seed
DO beams
Calculate IC and IT for each layer
Calcute and store flatwise and edewise MOE, thickness, density
END DO END DO
OUTPUT statistics:
beam thickness and density, flatwise and edgewise MOE,
regression parameters
(Program LVLBend continued)
DOUBLE PRECISION, DIMENSION(:,:), ALLOCATABLE :: Ee,Ef,Th,Dy EXTERNAL RNSET,RNNOF,RNUNF
INTEGER Iseed, I, J, N, Q, P, M INTEGER, DIMENSION(15) :: S
M=20 ! Number of Monte-Carlo Simulations N=20 ! Number of beams per simulation
! The names of the files to contain regressional info. and statistics for
! edgewise and flatwise MOE, thickness and density of the simulated beams FReg="C:\simres\LVLbend_Reg.txt"
FEedge="C:\simres\LVLbend_EDGE.txt"
FEflat="C:\simres\LVLbend_FLAT.txt"
FThick="C:\simres\LVLbend_THICK.txt"
FDens="C:\simres\LVLbend_DENS.txt"
W=5 ! Beam width (does not effect simulation results.) S = (/ 3,3,3,3,3,3,3,3,3,3,3,3,3,3,3 /) ! The layup
ALLOCATE ( Ee(M,N),Ef(M,N),Th(M,N),Dy(M,N) )
OPEN (UNIT=100,STATUS="OLD",file="c:\My documents\fortran\
&ranstream.txt")
! ranstream.txt contains a streem of random integers between 0 and 10000.
DO Q=1,M ! Executes M number of simulations WRITE (*,fmt="(I2,T13,\)"),Q
DO J=1,N ! Generates N beams for each simulation
! Reading and initializing the random seed READ (100,FMT="(E20.15)") A
Iseed=A
CALL RNSET(ISEED)
WRITE (*,fmt="(I3,\)"),J
DO I=1,10 ! This cycle is required because the NN = Sim("N",A,A) ! first simulated value is biased in
END DO ! the negative tail
! Generating layer properties DO I=1,15
Tho(I)=Orig(S(I)) Wei(I)=Tho(I)*Dty(S(I))
SELECT CASE (I) CASE (1,2,14,15)
(Program LVLBend contiunued) CASE (3:13)
Thi(I)=Core(Tho(I)) END SELECT
Ec(I)= CompHank(0.,0.,S(I)) * DENS( Tho(I),Thi(I),S(I) ) Et(I)= TenMOE(0.,S(I)) * DENS( Tho(I),Thi(I),S(I) ) END DO
! Moment of inertia in flatwise bending MP=SUM(THI)/2
! Calculating the flatwise and edgewise effective MOE Efl = SUM(Ec*Icf+Et*Itf) / Ifl
Eed = SUM(Ec*Thi+Et*Thi) / (2*SUM(Thi))
! Storing the flatwise and edgewise MOE, beam thickness and density Ee(Q,J)=Eed
Ef(Q,J)=Efl
(Program LVLBend continued)
END DO ! Ends 'NM' CLOSE (100)
PRINT *,""
! Saving the summary statistics of every simulation run;
! Separate files are created for each property.
OPEN (UNIT=200,STATUS="REPLACE",FILE=FReg) OPEN (UNIT=201,STATUS="REPLACE",FILE=FEedge) OPEN (UNIT=202,STATUS="REPLACE",FILE=FEflat) OPEN (UNIT=203,STATUS="REPLACE",FILE=FThick) OPEN (UNIT=204,STATUS="REPLACE",FILE=FDens) DO Q=1,M
Regpar=Reg( Ef(Q,:),Ee(Q,:) )
WRITE (200,FMT="(I3,3(F10.4))"),Q,Regpar(:) Stat = Sta(Ee(Q,:))
WRITE (201,FMT="(I3,5(F12.8))"),Q,Stat(:) Stat = Sta(Ef(Q,:))
WRITE (202,FMT="(I3,5(F12.8))"),Q,Stat(:) Stat = Sta(Th(Q,:))
WRITE (203,FMT="(I3,5(F12.8))"),Q,Stat(:) Stat = Sta(Dy(Q,:))
WRITE (204,FMT="(I3,5(F12.5))"),Q,Stat(:) END DO
CLOSE (200) CLOSE (201) CLOSE (202) CLOSE (203) CLOSE (204)
END PROGRAM LVLBend
The PSLBend program
Figure E.2 shows the flowchart for the program PSLBend, which simulates the bending MOE of PSL. This program is fairly similar to LVLBend, but is somewhat more complex, includes more random variables. Most notably, the calculation of the Ii values is more complicated, involving numerical integration. Instead of a layup, a species mix is specified. The first number in the variable Smix is the proportion of the species denoted by species code 1, the
START
Initialization of Species mix
Distribution parameters Cross-sectional dimensions
OUTPUT statistics:
ΣAi, density, MOE and ΣIi
flatwise and edgewise, regression parameters
END DO END DO
Calculate IC and IT for each strand
Calcute and store MOE and ΣIi
flatwise and edewise, ΣAi and density Generate horizontal and vertical
strand order
Generate horizontal and vertical strand position DO simulations
Generate for each strand:
Strand angle and orientation, Original and final thickness, Projected strand width, Species,Weight, ET and EC. (including the
densification effect) INPUT random seed
DO beams
Generate number of strands
ΣAi > A or ΣAi < 0.9485A
True
False
Table E.1 – The species-specific parameters specified in the Species module
Name Description Am Mean strand angle
As Standard deviation of strand angle Fm Average number of flakes per in2
Fs Standard deviation of number of flakes per in2 Om Mean strand deviation
Os Standard deviation of strand deviation Tm Mean final strand thickness (in)
Ts Standard deviation of final strand thickness (in)
second is that coded 2, etc. (e.g. for a pure yellow-poplar beam, Smix can contain the numbers {0, 0, 100}, since the code of yellow-poplar is 3.)
The parameters of the normal distributions for simulating the geometric parameters of the constituents are specified in the variables shown in Table E.1. R is the resolution of the numerical integration – increasing it makes the simulation more accurate but longer to run.
Program PSLBend USE Species USE Ifaces
CHARACTER (LEN=50) :: FReg,FEflat,FEedge,FIedge,FIflat,FArea,Fdens DOUBLE PRECISION Fm,Fs,Om,Os,Am,As,Thm,Ths,X,W,T
DOUBLE PRECISION L0,H,L,A,NN,Treal,Theta,Tot,f,fii,the DOUBLE PRECISION, DIMENSION(2) :: Is
DOUBLE PRECISION, DIMENSION(3) :: Regpar, Smix DOUBLE PRECISION, DIMENSION(5) :: Stat
DOUBLE PRECISION, DIMENSION(0:250) :: Ori,Ang,Thi,Ver,Wei,Rn,Hor,
& Ice,Ite,Icf,Itf,Eco,Ete,Tho,Lnt
DOUBLE PRECISION, DIMENSION(:,:), ALLOCATABLE :: Ee,Ef,Ied,Ifl,Ar,Dy INTEGER Iseed, I, J, NoS, N, Q, P, M, Temp, R
INTEGER, DIMENSION(0:250) :: Ho,Sps
EXTERNAL RNSET,RNNOF,RNUNF
! The names of the files to contain regressional info. and statistics for
! edgewise and flatwise MOE and moment of inertia,
! combined cross-section and density of the simulated beams
(Program PSLBend continued) FReg="C:\simres\PSLbend_Reg.txt"
FEedge="C:\simres\PSLbend_Eedge.txt"
FEflat="C:\simres\PSLbend_Eflat.txt"
FIedge="C:\simres\PSLbend_Iedge.txt"
FIflat="C:\simres\PSLbend_Iflat.txt"
FArea="C:\simres\PSLbend_Area.txt"
FDens="C:\simres\PSLbend_Dens.txt"
M=20 ! Number of Monte-Carlo Simulations N=20 ! Number of beams per simulation Smix=(/0,0,100/)
W=3 H=5.5 Fm=11.64 Fs=0.356 Om=2.991 Os=14.709 Am=0.03 As=4.668 Thm=0.08367 Ths=0.01542 R=10000
ALLOCATE (Ee(M,N),Ef(M,N),Ied(M,N),Ifl(M,N),Ar(M,N),Dy(M,N))
OPEN (UNIT=100,STATUS="OLD",FILE="c:\My documents\fortran\
&ranstream.txt")
! ranstream.txt contains a streem of random integers between 0 and 10000.
DO Q=1,M ! Executes M number of simulations WRITE (*,fmt="(I2,T13,\)"),Q
DO P=1,N ! Simulates N beams for each simulation
! Reading and initializing the random seed READ (100,FMT="(E20.15)") A
Iseed=A
CALL RNSET(ISEED)
WRITE (*,fmt="(I3,\)"),P
(Program PSLBend continued)
1 NN = W*H*Sim("N",Fm,Fs) ! Simulates the number of strands
NoS = INT(NN) ! Integer part (full width strands) L0 = NN-NoS ! Fractional part (width of strand nr 0.)
! Simulation of strand parameters DO I=0,NoS
Ori(I) = Sim( "N",Om,Os ) Ang(I) = Sim( "N",Am,As ) Sps(I) = SChoose(Smix)
Lnt(I) = 1/SQRT(COSD(Ang(I))**2*COSD(Ori(I))**2+SIND(Ori(I))**2) Tho(I) = Orig(Sps(I))`
5 Thi(I) = Sim("N",Thm,Ths)
Wei(I) = Tho(I)*Lnt(I)*Dty(Sps(I))
Treal=THI(I)*SQRT(SIND(Ori(I))**2*COSD(Ang(I))**2+COSD(Ori(I))**2) IF (Treal>Tho(I)) THEN ! This makes sure that the final thickness is
GO TO 5 ! Smaller than the original.
END IF
Lnt(0)=Lnt(0)*L0 ! The length and weight of the first strand Wei(0)=Wei(0)*L0 ! Is decreased according to L0
Tot=SUM(Thi(0:NoS)*Lnt(0:NoS)) ! The combined strand cross-section
! Beam rejections based on the combined cross-section IF (Tot>H*W) THEN
GO TO 1
ELSE IF (Tot<0.95*H*W) THEN GO TO 1
END IF
! Simulation of the horizontal order (vertical order is unchanged) DO I=0,NoS
(Program PSLBend continued)
! Simulation of the horizontal and vertical position
Ver(0) = (Thi(0)*Lnt(0)/Tot)*H/2-H/2
Hor(HO(0)) = (Thi(HO(0)*Lnt(HO(0)))/Tot)*W/2-W/2
DO I=1,NoS
Ver(I) = Ver(I-1) + (Thi(I)*Lnt(I)/Tot)*H
Hor(HO(I)) = Hor(HO(I-1)) + (Thi(HO(I))*Lnt(HO(I))/Tot)*W END DO
! Calculation of moment of inertia values horizontally and vertically DO I=0,NoS
Is = Numint(Ori(I),Ver(I),Thi(I),Lnt(I),H,R) Ice(I) = Is(1)
Ite(I) = Is(2)
Is = Numint(Ori(I)-90,Hor(I),Thi(I),Lnt(I),W,R) Icf(I) = Is(1)
Itf(I) = Is(2) END DO
! Calculating and storing flatwise and edgewise MOE and combined
! moment of inertia,the combined cross-section of the strands
! and beam density
Ee(Q,P)=SUM(Ice(0:NoS)*Eco(0:NoS)+Ite(0:NoS)*Ete(0:NoS))/(H**3*W/12) Ef(Q,P)=SUM(Icf(0:NoS)*Eco(0:NoS)+Itf(0:NoS)*Ete(0:NoS))/(H*W**3/12) Ied(Q,P)=SUM(Ice(0:NoS)+Ite(0:NoS))
Ifl(Q,P)=SUM(Icf(0:NoS)+Itf(0:NoS)) Ar(Q,P)=Tot
Dy(Q,P)=SUM(Wei(0:NoS))/(W*H)
END DO ! Ends 'N'
print *,""
END DO ! Ends 'M' CLOSE (100)
print *,""
! Saving the summary statistics of every simulation run;
! Separate files are created for each property.
(Program PSLBend continued)
DO Q=1,M
Regpar=Reg( Ef(Q,:),Ee(Q,:) )
WRITE (200,FMT="(I3,3(F10.4))"),Q,Regpar(:) Stat = Sta(Ee(Q,:))
WRITE (201,FMT="(I3,5(F12.8))"),Q,Stat(:) Stat = Sta(Ef(Q,:))
WRITE (202,FMT="(I3,5(F12.8))"),Q,Stat(:) Stat = Sta(Ied(Q,:))
WRITE (203,FMT="(I3,5(F12.8))"),Q,Stat(:) Stat = Sta(Ifl(Q,:))
WRITE (204,FMT="(I3,5(F12.8))"),Q,Stat(:) Stat = Sta(Ar(Q,:))
WRITE (205,FMT="(I3,5(F12.8))"),Q,Stat(:) Stat = Sta(Dy(Q,:))
WRITE (206,FMT="(I3,5(F12.5))"),Q,Stat(:) END DO
CLOSE (200) CLOSE (201) CLOSE (202) CLOSE (203) CLOSE (204)
end program PSLBend
The LVLComp program
Figure E.3 shows the flowchart for the program LVLComp, which simulates the orthotropic compression MOE of LVL in a certain direction. The program is similar to, though simpler than LVLBend. The variables Fip and Thp contain the ϕ’ and θ’ values associated with the compression load, respectively. The program outputs only the statistics for the simulated compression MOE.
Program LVLComp USE Species USE Ifaces
CHARACTER (LEN=50) :: FEcomp DOUBLE PRECISION X,T,A,Fip,Thp
DOUBLE PRECISION, DIMENSION(5) :: Stat
DOUBLE PRECISION, DIMENSION(:,:), ALLOCATABLE :: E EXTERNAL RNSET,RNNOF,RNUNF
INTEGER Iseed, I, J, N, Q, P, M, NN INTEGER, DIMENSION(15) :: S
M=20 ! Number of Monte Carlo Simulations N=10 ! Number of specimens per simulation
! The name of the files to contain regressional MOE statistics FEcomp="C:\simres\LVLcomp_G90_R0.txt"
Fip=90 Thp=0
S = (/ 3,3,3,3,3,3,3,3,3,3,3,3,3,3,3 /) ! The layup ALLOCATE ( E(M,N) )
OPEN (UNIT=100,STATUS="OLD",file="c:\My documents\fortran\
&ranstream.txt")
DO Q=1,M ! Executes M number of simulations WRITE (*,fmt="(I2,T13,\)"),Q
DO J=1,N ! Generates N specimens for each simulation
! Reading and initializing the random seed READ (100,FMT="(E20.15)") A
Iseed=A
CALL RNSET(Iseed)
WRITE (*,fmt="(I3,\)"),J
DO I=1,10 ! This cycle is required because the NN = Sim("N",A,A) ! first simulated value is biased in
END DO ! the negative tail
! Generating layer properties DO I=1,15
Tho(I)=Orig(S(I)) SELECT CASE (I) CASE (1,2,14,15)
Thi(I)=Face(Tho(I)) CASE (3:13)
Thi(I)=Core(Tho(I))
START
Calcute and store the beam’s compression MOE
END DO END DO
OUTPUT compression MOE statistics
STOP DO simulations
Generate for each layer:
Original thickness Layer thickness Eϕ’θ’ (including the
densification effect) INPUT random seed
DO beams Specification of the
layup, ϕ’ and θ’
Figure E.3 – Flowchart of the program LVLComp
(Program LVLComp continued)
! Calculating beam thickness and compression MOE T=SUM(THI)
E(Q,J)=SUM(MOE*THI)/T
END DO ! Ends 'N'
print *,""
END DO ! Ends 'NM' CLOSE (100)
print *,""
(Program LVLComp continued)
OPEN (UNIT=200,STATUS="REPLACE",FILE=FEcomp) DO Q=1,M
Stat = Sta(E(Q,:))
WRITE (200,FMT="(I3,5(F12.8))"),Q,Stat(:) END DO
CLOSE (200)
end program LVLComp
The PSLComp program
Figure E.4 shows the flowchart for the program PSLComp, which simulates the orthotropic compression MOE of PSL in a given direction. See some remarks at PSLbend and LVLComp.
Program PSLComp USE Species USE Ifaces
CHARACTER (LEN=50) :: FEcomp