Code Comments
Programming Forum and web based access to our favorite programming groups.Heres a fortran code to solve the equation of motion of very small particle in a low reynolds number flow. It is not working in all compilers, I have a working version of this code, but it doesnt work for all inputs. Can someone tell me what is wrong with this and fix it so that it works for a stokes number of 0.1 through 2. Stokes number is calculated in the program. You need to work backwards by playing with the input file and see for what values of input the stokes number is 0.1 through 2. my particle hardly moves. I have no way of even telling if my values are correct. What should I do? Pls help. my input file. it contains values of t_max, alpha, Rho, Rp, Rf, Uf in that order exactly example working input file values are: 600,0.05,2165,3.00E-07,5.20E-06,0.003 in that order. source code: IMPLICIT NONE REAL*8 DERIV_FX, DERIV_FY,N,RP,RF,RHO,UF,alpha,K,R0,U REAL*8 Nu,x,y,y0,THETA,DT,T0,UT0,PI,T,THETA0,St okes REAL*8 UR0,b,x1,VX,VY,R,UT_K,UR_K,THETA2,MP,R2, x_i INTEGER t_max, time REAL*8 n_t EXTERNAL FX, FY, DERIV_FX, DERIV_FY COMMON Nu, RP, RF, MP, UF, alpha, K, PI, Stokes, x_i CALL read_parameters (t_max, alpha, RHO, RP, RF, UF) PI = 3.1415 Nu = 0.0000181 MP = RHO * ((4.0/3.0)*PI*RP**3.0) Stokes = (2.0*RP) **2*RHO*UF/ (6*PI*Nu*2.0*RF) write (6,*) 'The Stokes is equal to:' write (6,*) Stokes b = RF * alpha**(-1.0/2.0) K = -0.5*log(alpha)-0.75+alpha-0.25*alpha**2.0 n_t = 1.0 5 y0 = n_t*b x= (-1.0)*(b**2.0 - y0**2.0)**(1.0/2.0) R0 = (x**2.0 + y0**2.0)**(1.0/2.0) THETA0 = atan(y0/x) THETA0= PI/2.0 + THETA0 DT = 0.000000001 T0 = 0.0 T =T0 THETA = THETA0 R = R0 Y= Y0 write (6,*) 'n_t=', n_t, 'y=', y U = UF UT0 = (U/(2.0*K))*(2.0*log(R/RF)+1+alpha-(RF**2.0/R**2.0) *(1.0-alpha/ 2.0)-(3.0*alpha/2.0)*(R**2.0/RF**2.0))*sin(THETA) UR0 = (U/(2.0*K))*(2.0*log(R/RF)-1.0+alpha+(RF**2.0/R**2.0) *(1.0- alpha/2.0)-(alpha/2.0)*(R**2.0/RF**2.0))*cos(THETA) VX= UT0 VY = UR0 OPEN (11, FILE='PARTICLEOUT.DAT') WRITE (11,*) 'OUTPUT VARIABLES = time,THETA, R, x, y, VT, UT,VR, UR' DO 7 time = 1, t_max x_i =x CALL KUTTA(x,y,VX,VY,T,DT,FX,FY) THETA2 = atan(y/x) R2 = (x**2.0 + y**2.0)**(1.0/2.0) UT_K = (U/(2.0*K))*(2.0*log(R2/RF)+1.0+alpha-(RF**2.0/R2**2.0)*(1.0- alpha/2.0)-3.0*alpha/2.0*(R2**2.0/RF**2.0))*sin(THETA2) UR_K = (U/(2.0*K))*(2.0*log(R2/RF)-1.0+alpha+(RF**2.0/R2**2.0)*(1.0- alpha/2.0)-alpha/2.0*(R2**2.0/RF**2.0))*cos(THETA2) WRITE(11,10) time,THETA2,R2,x,y,VX,UT_K,VY,UR_K 10 FORMAT(I1,9E15.6) if ((x**2+y**2) .LE. (RF+RP)**2) then write (6,*) 'The Stokes is equal to:' write (6,*) Stokes stop endif if (n_t .LE. 8.0E-8) then close (11) write (6,*) 'Neglegible single fiber efficiency under these conditions' stop endif if ((x-x_i) .LE. 0.0) then n_t = n_t - (n_t/200.0) close (11) goto 5 stop endif 7 CONTINUE CLOSE (11) END SUBROUTINE KUTTA(x,y,VX,VY,T,DT,FX,FY) IMPLICIT NONE REAL*8 D1X,D1Y,DT,VX,VY,D1U,D1V,D2X,D2Y,D2U,D2V REAL*8 D3X,D3Y,D3U,D3V,D4X,D4Y,D4U,D4V,T,x,y REAL*8 DERIV_FX, DERIV_FY DERIV_FX = 0.0 DERIV_FY =0.0 D1X = DT * VX D1Y = DT * VY CALL FX(x,y,VX,DERIV_FX) CALL FY(x,y,VY,DERIV_FY) D1U = DT * DERIV_FX D1V= DT * DERIV_FY D2X = DT * (VX+D1U/2.0) D2Y = DT * (VY+D1V/2.0) CALL FX(x+D1X/2.0, y+D1Y/2.0, VX+D1U/ 2.0, DERIV_FX) CALL FY(x+D1X/2.0, y+D1Y/2.0, VY+D1V/ 2.0, DERIV_FY) D2U = DT * DERIV_FX D2V = DT * DERIV_FY D3X = DT * (VX+D2U/2.0) D3Y = DT * (VY+D2V/2.0) CALL FX(x+D2X/2.0, y+D2Y/2.0, VX+D2U/ 2.0, DERIV_FX) CALL FY(x+D2X/2.0, y+D2Y/2.0, VY+D2V/ 2.0, DERIV_FY) D3U = DT * DERIV_FX D3V = DT * DERIV_FY D4X = DT * (VX+D3U) D4Y = DT * (VY+D3V) CALL FX(x+D3X, y+D3Y, Vx+D3U, DERIV_FX) CALL FY(x+D3X, y+D3Y, VY+D3V, DERIV_FY) D4U = DT * DERIV_FX D4V= DT * DERIV_FY T = T + DT X= x + (D1X + 2.0*D2X + 2.0*D3X + D4X) / 6.0 y =y + (D1Y + 2.0*D2Y + 2.0*D3Y + D4Y) / 6.0 VX = VX + (D1U + 2.0*D2U + 2.0*D3U + D4U) / 6.0 VY = VY + (D1V + 2.0*D2V + 2.0*D3V + D4V) / 6.0 RETURN END SUBROUTINE FX (x, y, VX, DERIV_FX) IMPLICIT NONE REAL*8 UT_K,Sk,K,Nu,alpha,y,RP,RF,MP,UF,x,VX REAL*8 Stokes,DERIV_FX,PI,THETA,R,x_i,U COMMON Nu,RP,RF,MP,UF,alpha,K,PI,Stokes,x_i R= (x**2.0 + y**2.0)**(1.0/2.0) U = UF THETA = atan(y/x) UT_K = (U/(2.0*K))*(2.0*log(R/RF) +1.0+alpha-(RF**2.0/R**2.0)*(1.0-alpha/2.0)-(3.0*alpha/2.0)*(R**2.0/ RF**2.0))*sin(THETA) Sk = 6.0*PI*Nu*RP/MP DERIV_FX = Sk*(UT_K-VX) RETURN END SUBROUTINE FY(x, y, VY, DERIV_FY) IMPLICIT NONE REAL*8 UR_K, Sk, K, Nu, alpha, y, RP, RF, MP, UF, THETA, VY REAL*8 Stokes,DERIV_FY,PI,x,R,x_i,U COMMON Nu,RP,RF,MP,UF,alpha,K,PI,Stokes,x_i R= (x**2.0 + y**2.0)**(1.0/2.0) U = UF THETA = atan(y/x) UR_K=( (U/(2.0*K))*(2.0*log(R/RF)-1.0+alpha+ (RF**2.0/R**2.0)*(1.0-alpha/2.0)-(alpha/2.0)*(R**2.0/ RF**2.0))*cos(THETA)) Sk = 6.0*PI*Nu*RP/MP DERIV_FY = Sk* (UR_K-VY) RETURN END SUBROUTINE read_parameters (t_max,alpha, RHO, RP, RF, UF) implicit none integer t_max real*8 alpha,RHO,RP,RF,UF open (10, file='particle.par') read(10,*) t_max read(10,*) alpha read(10,*) RHO read(10,*) RP read(10,*) RF read(10,*) UF close (10) write (6,*) '*** Parameters read from file particle.par' write (6,*) '***' return end
Post Follow-up to this messageIn article <2536c6bf-4017-43c3-aa00-5d9eb8dd305e@u72g2000hsf.googlegroups.co m>, madhu0319@gmail.com writes: > Heres a fortran code to solve .... > source code: > > IMPLICIT NONE > REAL*8 DERIV_FX, > DERIV_FY,N,RP,RF,RHO,UF,alpha,K,R0,U > REAL*8 > Nu,x,y,y0,THETA,DT,T0,UT0,PI,T,THETA0,St okes > REAL*8 > UR0,b,x1,VX,VY,R,UT_K,UR_K,THETA2,MP,R2, x_i > INTEGER t_max, time > REAL*8 n_t > EXTERNAL FX, FY, DERIV_FX, DERIV_FY > COMMON Nu, RP, RF, MP, UF, alpha, K, > PI, Stokes, x_i > CALL read_parameters (t_max, alpha, > RHO, RP, RF, UF) > PI = 3.1415 > Nu = 0.0000181 The problems in this code are too numerous to expend the effort until you actual post code that has a just to compile with any compiler. -- steve
Post Follow-up to this messagemadhu0319@gmail.com wrote: > Heres a fortran code to solve the equation of motion of very small > particle in a low reynolds number flow. It is not working in all > compilers, I have a working version of this code, but it doesnt work > for all inputs. Can someone tell me what is wrong with this and fix it > so that it works for a stokes number of 0.1 through 2. Stokes number > is calculated in the program. You need to work backwards by playing > with the input file and see for what values of input the stokes number > is 0.1 through 2. > > my particle hardly moves. I have no way of even telling if my values > are correct. What should I do? Pls help. > my input file. it contains values of t_max, alpha, Rho, Rp, Rf, Uf in > that order exactly > example working input file values are: > 600,0.05,2165,3.00E-07,5.20E-06,0.003 in that order. > > source code: > > IMPLICIT NONE > REAL*8 DERIV_FX, > DERIV_FY,N,RP,RF,RHO,UF,alpha,K,R0,U > REAL*8 > Nu,x,y,y0,THETA,DT,T0,UT0,PI,T,THETA0,St okes > REAL*8 > UR0,b,x1,VX,VY,R,UT_K,UR_K,THETA2,MP,R2, x_i > INTEGER t_max, time > REAL*8 n_t > EXTERNAL FX, FY, DERIV_FX, DERIV_FY > COMMON Nu, RP, RF, MP, UF, alpha, K, > PI, Stokes, x_i > CALL read_parameters (t_max, alpha, > RHO, RP, RF, UF) > PI = 3.1415 I don't know enough about your problem to debug it. But, here's a simple problem. You've declared all of you variables to be real*8, which is usually called double precision. But, you specify all of your constants as single precision. 3.14159 is a single precision value. The fact that it is assigned to a double precision variable does not make it more precise; values are usually converted from single to double by appending zeros in the low order bits. It doesn't make much difference for PI, since you only specify a few digits. For any number that isn't an integer or a "power of 2" fraction like .5 or .25, specify the value in double precision. The way to do that, consistent with your coding style, is to write 3.1415D0 Better yet to look up a few more digits for PI, it can't hurt to be more accurate ;) . > Nu = 0.0000181 > MP = RHO * ((4.0/3.0)*PI*RP**3.0) This is a potentially worse problem. The expression (4.0/3.0) will be computed as approximately 1.333333 and then converted to double precision and become approximately 1.3333330000000. You probably wanted 1.3333333333333. For expressions like (4.0/3.0), make sure at least one of the values is double precision. Best to do (4.0D0/3.0D0) . Secondly, you almost never want to use expressions like RP**3.0 where the exponent is a real number. Most compilers will do that with log and exp functions and lose accuracy. Better to do exponentiation as RP**3 when the exponent is an integer. > Stokes = (2.0*RP) **2*RHO*UF/ > (6*PI*Nu*2.0*RF) > write (6,*) 'The Stokes is equal to:' > write (6,*) Stokes > b = RF * alpha**(-1.0/2.0) Personally, I'd use the SQRT function, rather than exponentiation to .5. It's likely to be faster, more accurate, and much easier for humans to read. Dick Hendrickson > K = > -0.5*log(alpha)-0.75+alpha-0.25*alpha**2.0 > n_t = 1.0 > 5 y0 = n_t*b > x= (-1.0)*(b**2.0 - > y0**2.0)**(1.0/2.0) > R0 = (x**2.0 + y0**2.0)**(1.0/2.0) > THETA0 = atan(y0/x) > > > THETA0= PI/2.0 + THETA0 > DT = 0.000000001 > T0 = 0.0 > T =T0 > THETA = THETA0 > R = R0 > Y= Y0 > write (6,*) 'n_t=', n_t, 'y=', y > U = UF > > UT0 = (U/(2.0*K))*(2.0*log(R/RF)+1+alpha-(RF**2.0/R**2.0) *(1.0-alpha/ > 2.0)-(3.0*alpha/2.0)*(R**2.0/RF**2.0))*sin(THETA) > UR0 = (U/(2.0*K))*(2.0*log(R/RF)-1.0+alpha+(RF**2.0/R**2.0) *(1.0- > alpha/2.0)-(alpha/2.0)*(R**2.0/RF**2.0))*cos(THETA) > VX= UT0 > VY = UR0 > OPEN (11, FILE='PARTICLEOUT.DAT') > WRITE (11,*) 'OUTPUT VARIABLES = > time,THETA, R, x, y, VT, UT,VR, UR' > DO 7 time = 1, t_max > x_i =x > CALL KUTTA(x,y,VX,VY,T,DT,FX,FY) > THETA2 = atan(y/x) > R2 = (x**2.0 + y**2.0)**(1.0/2.0) > UT_K = (U/(2.0*K))*(2.0*log(R2/RF)+1.0+alpha-(RF**2.0/R2**2.0)*(1.0- > alpha/2.0)-3.0*alpha/2.0*(R2**2.0/RF**2.0))*sin(THETA2) > UR_K = (U/(2.0*K))*(2.0*log(R2/RF)-1.0+alpha+(RF**2.0/R2**2.0)*(1.0- > alpha/2.0)-alpha/2.0*(R2**2.0/RF**2.0))*cos(THETA2) > WRITE(11,10) > time,THETA2,R2,x,y,VX,UT_K,VY,UR_K > > 10 FORMAT(I1,9E15.6) > > if ((x**2+y**2) .LE. (RF+RP)**2) then > write (6,*) 'The Stokes is equal to:' > write (6,*) Stokes > stop > endif > > if (n_t .LE. 8.0E-8) then > close (11) > write (6,*) 'Neglegible single fiber > efficiency under these conditions' > > stop > endif > > if ((x-x_i) .LE. 0.0) then > n_t = n_t - (n_t/200.0) > close (11) > goto 5 > stop > endif > 7 CONTINUE > CLOSE (11) > END > > SUBROUTINE KUTTA(x,y,VX,VY,T,DT,FX,FY) > IMPLICIT NONE > REAL*8 > D1X,D1Y,DT,VX,VY,D1U,D1V,D2X,D2Y,D2U,D2V > REAL*8 > D3X,D3Y,D3U,D3V,D4X,D4Y,D4U,D4V,T,x,y > REAL*8 DERIV_FX, DERIV_FY > > DERIV_FX = 0.0 > DERIV_FY =0.0 > > D1X = DT * VX > D1Y = DT * VY > > CALL FX(x,y,VX,DERIV_FX) > CALL FY(x,y,VY,DERIV_FY) > > D1U = DT * DERIV_FX > D1V= DT * DERIV_FY > > D2X = DT * (VX+D1U/2.0) > D2Y = DT * (VY+D1V/2.0) > > CALL FX(x+D1X/2.0, y+D1Y/2.0, VX+D1U/ > 2.0, DERIV_FX) > CALL FY(x+D1X/2.0, y+D1Y/2.0, VY+D1V/ > 2.0, DERIV_FY) > > D2U = DT * DERIV_FX > D2V = DT * DERIV_FY > > D3X = DT * (VX+D2U/2.0) > D3Y = DT * (VY+D2V/2.0) > > CALL FX(x+D2X/2.0, y+D2Y/2.0, VX+D2U/ > 2.0, DERIV_FX) > CALL FY(x+D2X/2.0, y+D2Y/2.0, VY+D2V/ > 2.0, DERIV_FY) > > D3U = DT * DERIV_FX > D3V = DT * DERIV_FY > > D4X = DT * (VX+D3U) > D4Y = DT * (VY+D3V) > > CALL FX(x+D3X, y+D3Y, Vx+D3U, > DERIV_FX) > CALL FY(x+D3X, y+D3Y, VY+D3V, > DERIV_FY) > > D4U = DT * DERIV_FX > D4V= DT * DERIV_FY > > T = T + DT > X= x + (D1X + 2.0*D2X + 2.0*D3X + > D4X) / 6.0 > y =y + (D1Y + 2.0*D2Y + 2.0*D3Y + > D4Y) / 6.0 > VX = VX + (D1U + 2.0*D2U + 2.0*D3U + > D4U) / 6.0 > VY = VY + (D1V + 2.0*D2V + 2.0*D3V + > D4V) / 6.0 > > RETURN > END > > SUBROUTINE FX (x, y, VX, DERIV_FX) > > IMPLICIT NONE > REAL*8 > UT_K,Sk,K,Nu,alpha,y,RP,RF,MP,UF,x,VX > REAL*8 > Stokes,DERIV_FX,PI,THETA,R,x_i,U > COMMON > Nu,RP,RF,MP,UF,alpha,K,PI,Stokes,x_i > > R= (x**2.0 + y**2.0)**(1.0/2.0) > U = UF > > THETA = atan(y/x) > UT_K = (U/(2.0*K))*(2.0*log(R/RF) > +1.0+alpha-(RF**2.0/R**2.0)*(1.0-alpha/2.0)-(3.0*alpha/2.0)*(R**2.0/ > RF**2.0))*sin(THETA) > Sk = 6.0*PI*Nu*RP/MP > > DERIV_FX = Sk*(UT_K-VX) > > RETURN > END > > SUBROUTINE FY(x, y, VY, DERIV_FY) > > IMPLICIT NONE > REAL*8 UR_K, Sk, K, Nu, alpha, y, RP, RF, MP, > UF, THETA, VY > REAL*8 Stokes,DERIV_FY,PI,x,R,x_i,U > COMMON Nu,RP,RF,MP,UF,alpha,K,PI,Stokes,x_i > > R= (x**2.0 + y**2.0)**(1.0/2.0) > U = UF > > > THETA = atan(y/x) > UR_K=( (U/(2.0*K))*(2.0*log(R/RF)-1.0+alpha+ > (RF**2.0/R**2.0)*(1.0-alpha/2.0)-(alpha/2.0)*(R**2.0/ > RF**2.0))*cos(THETA)) > > Sk = 6.0*PI*Nu*RP/MP > > > DERIV_FY = Sk* (UR_K-VY) > > RETURN > END > > SUBROUTINE read_parameters (t_max,alpha, RHO, > RP, RF, UF) > implicit none > integer t_max > real*8 alpha,RHO,RP,RF,UF > open (10, file='particle.par') > > read(10,*) t_max > > read(10,*) alpha > > read(10,*) RHO > > read(10,*) RP > > read(10,*) RF > > read(10,*) UF > > close (10) > > write (6,*) '*** Parameters read from file > particle.par' > write (6,*) '***' > > return > end > >
Post Follow-up to this messageOn Mar 26, 8:57=A0pm, madhu0...@gmail.com wrote: > Heres a fortran code to solve the equation of =A0motion of very small > particle in a low reynolds number flow. It is not working in all > compilers, I have a working version of this code, but it doesnt work > for all inputs. Can someone tell me what is wrong with this and fix it > so that it works for a stokes number of 0.1 through 2. Stokes number > is calculated in the program. You need to work backwards by playing > with the input file and see for what values of input the stokes number > is 0.1 through 2. > > my particle hardly moves. I have no way of even telling if my values > are correct. What should I do? Pls help. > my input file. it contains values of t_max, alpha, Rho, Rp, Rf, Uf in > that order exactly > example working input file values are: > 600,0.05,2165,3.00E-07,5.20E-06,0.003 in that order. > > source code: > > IMPLICIT NONE > =A0 =A0 =A0 =A0 =A0 =A0 =A0 =A0 =A0 =A0 =A0 =A0 =A0 =A0 =A0 =A0 REAL*8 DER=[/color ] IV_FX, > DERIV_FY,N,RP,RF,RHO,UF,alpha,K,R0,U > =A0 =A0 =A0 =A0 =A0 =A0 =A0 =A0 =A0 =A0 =A0 =A0 =A0 =A0 =A0 =A0 REAL*8 > Nu,x,y,y0,THETA,DT,T0,UT0,PI,T,THETA0,St okes > =A0 =A0 =A0 =A0 =A0 =A0 =A0 =A0 =A0 =A0 =A0 =A0 =A0 =A0 =A0 =A0 REAL*8 > UR0,b,x1,VX,VY,R,UT_K,UR_K,THETA2,MP,R2, x_i > =A0 =A0 =A0 =A0 =A0 =A0 =A0 =A0 =A0 =A0 =A0 =A0 =A0 =A0 =A0 =A0 INTEGER t_=[/color ] max, time [snip rest of listing] 1. Is this fixed format or free format source code? Your compiler may truncate lines beyond 72 or 132 characters without warning. 2. Have you used an external diagnostic program to examine your code like FTNCHEK? 3. Have you used your compiler's diagnostic options? - e
Post Follow-up to this messageOn Mar 26, 8:19 pm, ka...@troutmask.apl.washington.edu (Steven G. Kargl) wrote: > In article <2536c6bf-4017-43c3-aa00-5d9eb8dd3...@u72g2000hsf.googlegroups. com>, > madhu0...@gmail.com writes: > > > > .... > > > The problems in this code are too numerous to expend > the effort until you actual post code that has a just > to compile with any compiler. > > -- > steve Sorry, I did not really understand what u said. I do have an executable file for this code that works, but i dont know how to post it here..What do I do abt this program, is it fixable at all?
Post Follow-up to this messageThankyou, is there any way i can attach files here? I have an exe file for the code above, but i dont know if it returns correct values, i dont anything abt coding. i am a zero. this is another guys work that i have to improve. i have to get better results and submit it..is this program feasible in a time frame of 20days?or is it entirely wrong. heres what is going on in the program.. What basically we are doing is tracking a particles trajectory and checking if it=A2s hitting a fiber filter that we have. It=A2s a low Reynolds number flow. That program is not working at very small stokes number and at every particle velocity. So I have to alter the program. For the particle trajectory we are solving an equation dv/dt=3D 18 * =E7*(U-V) / (2rp )2 =F1p ---------------------------------------- equation 1 where =E7 =3D (2*RP)^2*RHO*UF/ (6*PI*Nu*2*RF) This is also the equation 4.19 of the textbook that the program refers to while solving. The program goes like this 1. Reading Initialization parameters namely maximum number of iterations, packing density, particle density, particle radius, fiber radius and free stream fluid velocity. 2. Estimate the release point using following relationship: x=3D (b2- y02) ^0.5 and start iterations. Here y0 =3D ni*b where b is the radius of the imaginary fluid surrounding the fiber and ni is a random number designating the initial position of the particle and 0< ni < 1. B=3D RF*alpha^ (-1/2) 3. Start with ni=3D 1 and keep decreasing it by 1% till the following condition is satisfied 4. condition x^2 + y^2 <=3D (rf+rp)^2. Subscript f refers to fiber and Subscript p refers to particle. Here x and y are the coordinates obtained during the execution of the Runge Kutta 4 method for equation 1. 5. Release particle from the kuwabara cell. 6. Print output when the particle touches the fiber. For each iteration we solve equation 1 for X and Y coordinates points and the velocity components using the RK 4 method. Print x and y coordinates the velocity components after each iteration. 7. The coordinates are printed when the particle touches the fiber. 8. RK 4 method: This is used to calculate the value of of say a variable say y. By using its present value say yn, we can calculate the next value of y that is y (n+1). In our present program since we have to calculate x, y, VX, VY, T: We have 4 sets of k values by which I mean Our k1, k2, k3, k4 will be. So we have to solve that equation for x, y , VX. VY, . T is calculated by simple increment of time step that is T=3DT+DT. D1x D1y D1u D1v D2x D2y D2u D2v D3x D3y D3u D3v D4x D4y D4u D4v D1X actually means change in position in a certain time. Its Dx=3D Dt*VX ( change in time multiplied by the velocity X which means change in position ) . U can look at it this way : DX/DT =3D VX (distance in X direction over time is X velocity). Same is the case with D1Y. D1Y=3D Dt*VY. D1U is the change in velocity in the X direction which will be acceleration. We know that F=3DM*a. Which means acceleration a=3D change in velocity =3D F/M. Implies change in velocity in the Xdirection is (forces in Xdirection over M). So D1U/DT=3D FX or D1U=3D DT*FX. Since F is a sum of all forces F is defined as a function of x,y, u, v etc Calculate FX, Fy in separate functions to make things easier . Another reason FX and Fy are calculated separate is that they represent the RHS of equation of motion.
Post Follow-up to this messagei did not use any diagnostic program. i dont even know anything abt it. i dont know anything abt fortran or coding..this is my project at school thats due in 20days..cani run this code on any compiler, i mean is it compatible with any? shud it return the same values on any compiler? pls bear with my silly questions, i am very stressed out abt this..
Post Follow-up to this messageIn article <4e374a6d-1203-40a2-8b15-216aa141ea8d@x41g2000hsb.googlegroups.co m>, madhu0319@gmail.com writes: > On Mar 26, 8:19 pm, ka...@troutmask.apl.washington.edu (Steven G. > Kargl) wrote: > > Sorry, I did not really understand what u said. I do have an > executable file for this code that works, but i dont know how to post > it here..What do I do abt this program, is it fixable at all? I'm not interested in your executable. The code you posted isn't in a proper form, so no Fortran compiler can compile it. You need to post Fortran code that is in a *readable* and standard conforming form. I suggest that you get a few books from your local library on Fortran programming. You have 20 days to read the material, and perhaps learn something. -- Steve
Post Follow-up to this messageOn Mar 26, 10:08=A0pm, madhu0...@gmail.com wrote: > i did not use any diagnostic program. i dont even know anything abt > it. i dont know anything abt fortran or coding..this is my project at > school thats due in 20days..cani run this code on any compiler, i mean > is it compatible with any? shud it return the same values on any > compiler? > > pls bear with my silly questions, i am very stressed out abt this.. You CAN get your program to compile if you 1. remove leading spaces on each line 2. change EXTERNAL FX, FY, DERIV_FX, DERIV_FY to EXTERNAL FX, FY 3. compile as free format source DERIV_FX and DERIV_FY are variables NOT functions passed as arguments to routines so there is no need to make them EXTERNAL. Making the numerical changes suggested by others DOES not make any difference in the output except for trailing digits for the value of Stokes. When I compiled and ran the program it printed: *** Parameters read from file particle.par *** The Stokes is equal to: 0.0006589939662307433 n_t=3D 1. y=3D 0.000023255106965997813 n_t=3D 0.995 y=3D 0.000023138831431167825 n_t=3D 0.990025 y=3D 0.000023023137274011985 n_t=3D 0.985074875 y=3D 0.000022908021587641927 n_t=3D 0.9801495006250001 y=3D 0.000022793481479703716 n_t=3D 0.9752487531218751 y=3D 0.0000226795140723052 n_t=3D 0.9703725093562657 y=3D 0.000022566116501943672 n_t=3D 0.9655206468094844 y=3D 0.000022453285919433955 n_t=3D 0.960693043575437 y=3D 0.000022341019489836785 Note that your value of Stokes is NOT between 0.1 and 2.0 as desired. Maybe better input data will help. Maybe better values for other constants will help. Check the program listing against your formulas and equations. Verify that your program is correctly transcribed. Insert WRITE statements in the program to display the value of variables as the program runs. Compare these to manual calculations. What are reasonable ranges for the input variables? Try combinations of LOW, MID and HIGH values for the input variables and see what you get. At this point there is no substitute for pouring over your program listing and output with paper, pencil and a calculator. With MINOR changes the program DOES compile. You don't have obvious errors involving mismatched arguments to subroutines or common blocks that do not match. You may have other programming errors logical errors bad numerical behavior. Other than the "Non-Fortran" flavor of your program, there is nothing obviously wrong that jumps out except goto 5 followed by STOP. :-). Good luck. -- e
Post Follow-up to this messageIf the program becomes exhausted and can't make it to 5, maybe it will just fall into STOP. Good suggestions. Reminds me of my old old old old college days. On Wed, 26 Mar 2008 21:21:53 -0700 (PDT), e p chandler <epc8@juno.com> wrote: >Other than the "Non-Fortran" flavor of your program, there is nothing >obviously wrong that jumps out except > >goto 5 >followed by >STOP. :-).
Post Follow-up to this messagePowered by vBulletin
Copyright 2000-2006 Jelsoft Enterprises Limited.