{[ promptMessage ]}

Bookmark it

{[ promptMessage ]}

CE 311K Lab12

# CE 311K Lab12 - largest=abs(EIvpp Sheet1 Page 2 max_loc =...

This preview shows pages 1–2. Sign up to view the full content.

Sheet1 Page 1 ! CE311K Introduction to Computer Methods ! Lab12 - Differentiation Program Lab12 implicit none real x, x0, x1, x2, x3, x4, x5, x6,L,q,a,b,scale real h, smallest, max_loc, max_mom integer num, i ! Open output files Open (unit=6, file='dif1.dat', status='unknown') Open (unit=7, file='dif2.dat', status='unknown') ! Initialization of variables h=0.001 L=1. q=8000. smallest = 1000 ! Write headers for output files Write (7,60) Write (7,50) Write (6,60) Write (6,50) 50 Format(50('=')) 60 Format(3x,'Distance (m)', 6x, 'Bending Moment (N-m)') ! Calculation of bending moment using 1st-order differentiation call Differ1(L,h,q) call Differ2(L,h,q) End subroutine Differ1(L,h,q) implicit none real x, x0, x1, x2, x3, x4,x5,x6,h,max_mom,max_loc,largest,L,q,EIv,EIvpp,EIvppp integer i,num EIv(x)=q*L*x**3/12.-q*x**4/24.-q*L**3*x/24. num=int(L/h) largest=0. Do i=2, num-2 x0 = h*(i-2) x1 = x0+h x2 = x1+h x3 = x2+h x4 = x3+h EIvpp = (EIv(x3)-2*EIv(x2)+EIv(x1))/(h**2) EIvppp = (EIv(x4)-2*EIv(x3)+2*EIv(x1)-EIv(x0))/(2*h**3) If(abs(EIvpp).gt.largest) then

This preview has intentionally blurred sections. Sign up to view the full version.

View Full Document
This is the end of the preview. Sign up to access the rest of the document.

Unformatted text preview: largest=abs(EIvpp) Sheet1 Page 2 max_loc = x2 max_mom = largest Endif Write(6,5) x2, EIvpp Enddo Write(6,150) max_mom, max_loc,L 5 format(f5.4,5x,f8.3) 150 Format (2/,'Max bending moment is ',f8.2,' and occurs at ', f8.2,' for L=',f3.1) return end subroutine Differ2(L,h,q) implicit none real x, x0, x1, x2, x3, x4,x5,x6,h,max_mom,max_loc,largest,L,q,EIv,EIvpp,EIvppp integer i,num EIv(x)=q*L*x**3/12.-q*x**4/24.-q*L**3*x/24. num=int(L/h) largest=0. Do i=3, num-3 x0 = h*(i-3) x1 = x0+h x2 = x1+h x3 = x2+h x4 = x3+h x5 = x4+h x6 = x5+h EIvpp = (-EIv(x5)+16*EIv(x4)-30*EIv(x3)+16*EIv(x2)-EIv(x1))/(12*h**2) EIvppp = (-EIv(x6)+8*EIv(x5)-13*EIv(x4)+13*EIv(x2)-8*EIv(x1)+EIv(x0))/(8*h**3) If(abs(EIvpp).gt.largest) then largest=abs(EIvpp) max_loc = x3 max_mom = largest Endif Write(7,5) x3, EIvpp Enddo Write(7,150) max_mom, max_loc,L 5 format(f5.4,5x,f8.3) 150 Format (2/,'Max bending moment is ',f8.2,' and occurs at ', f8.2,' for L=',f3.1) return end...
View Full Document

{[ snackBarMessage ]}