 |
forums.silverfrost.com Welcome to the Silverfrost forums
|
View previous topic :: View next topic |
Author |
Message |
Srabon
Joined: 22 Feb 2011 Posts: 14
|
Posted: Sun Mar 04, 2012 8:54 pm Post subject: problem with coding pseudocode |
|
|
I am trying to code the following pseudo-code to approximate the solution to the 1D linear convection as in the figure:
http://imageshack.us/photo/my-images/708/fig1s.jpg/
Here is what I have done so far,
Code: | PROGRAM one_d_linear_convection
IMPLICIT NONE
INTEGER::i,it,nx,nt,k
DOUBLE PRECISION::dx,dt,c
DOUBLE PRECISION,DIMENSION(1:20)::u,un
nx=20
nt=50
dt=0.01
c=1.0
dx=2./(nx-1.)
!initial condition
DO i = 1,nx
IF (i*dx>=0.5 .and. i*dx<=1) THEN
u(i) = 2
ELSE
u(i) = 1
ENDIF
! WRITE(*,*)i,u(i)
ENDDO
!Finite Difference
DO it=1,nt
DO k=0,nx-1
un(k)=u(k)
DO i=2,nx-1
u(i) = un(i) - c*dt/dx*(un(i)-un(i-1))
WRITE(*,*)i,u(i)
ENDDO
ENDDO
ENDDO
END
|
My problem is I cant figure a way out to see the arrays and the u value at x=1 and 20. Any help would be highly appreciated.
Thanks! |
|
Back to top |
|
 |
JohnCampbell
Joined: 16 Feb 2006 Posts: 2615 Location: Sydney
|
Posted: Sun Mar 04, 2012 11:52 pm Post subject: |
|
|
Do these changes get closer to what you want Code: | PROGRAM one_d_linear_convection
IMPLICIT NONE
INTEGER::i,it,nx,nt,k
DOUBLE PRECISION::dx,dt,c, x
DOUBLE PRECISION,DIMENSION(1:20)::u,un
!
nx=20
nt=50
dt=0.01
c=1.0
dx=2./(nx-1.)
!initial condition
DO i = 1,nx
x = i*dx
IF (x>=0.5 .and. x<=1.0) THEN
u(i) = 2
ELSE
u(i) = 1
ENDIF
! WRITE(*,*)i,u(i)
END DO
!Finite Difference
DO it=1,nt
!
DO k=1,nx ! 0,nx-1
un(k)=u(k)
ENDDO
!
DO i=2,nx-1
u(i) = un(i) - c*dt/dx*(un(i)-un(i-1))
ENDDO
!
DO i=1,nx
WRITE(*,*) it,i,u(i)
ENDDO
!
ENDDO
END |
You can import the output from the program into excel and chart the results.
You now have flow in only one direction ?
why not have :
u(i) = un(i) - c*dt/dx*(un(i)-un(i-1)) &
- c*dt/dx*(un(i)-un(i+1)) |
|
Back to top |
|
 |
JohnCampbell
Joined: 16 Feb 2006 Posts: 2615 Location: Sydney
|
Posted: Mon Mar 05, 2012 5:13 am Post subject: |
|
|
A better solution for more itterations, printing fewer results and allowing flow both ways could be Code: | PROGRAM one_d_linear_convection
IMPLICIT NONE
INTEGER :: i,it,nx,nt,k
REAL*8 :: dx,dt,c, x
REAL*8,DIMENSION(40) :: u,un
!
nx = size (u)
dt = 0.001 ! use a smaller time step
nt = 501
c = 1.0
dx = 2./(nx-1.)
!
! initial condition
DO i = 1,nx
x = i*dx
IF (x>=0.5 .and. x<=1.0) THEN
u(i) = 2
ELSE
u(i) = 1
ENDIF
! WRITE(*,*)i,u(i)
END DO
!
! Finite Difference
DO it=1,nt
!
DO k=1,nx ! 0,nx-1
un(k)=u(k)
ENDDO
!
DO i=2,nx-1
u(i) = un(i) - c*dt/dx/2.0 * ( (un(i)-un(i-1)) + (un(i)-un(i+1)) )
ENDDO
!
if (mod(it,100)== 1) then ! only print every 0.1 seconds
DO i=1,nx
WRITE(*,*) it,i,u(i)
ENDDO
end if
!
ENDDO
!
END |
|
|
Back to top |
|
 |
Srabon
Joined: 22 Feb 2011 Posts: 14
|
Posted: Mon Mar 05, 2012 6:55 am Post subject: |
|
|
That was so nice of you, John.
Thank you so much for the codes! |
|
Back to top |
|
 |
Srabon
Joined: 22 Feb 2011 Posts: 14
|
Posted: Mon Apr 09, 2012 1:21 pm Post subject: |
|
|
Code: | PROGRAM odlc
IMPLICIT NONE
INTEGER::i,it
INTEGER,PARAMETER::nx=20,nt=50
DOUBLE PRECISION,PARAMETER::dx=2./(nx-1.),dt=0.01,c=1.0
DOUBLE PRECISION,DIMENSION(nt,nx)::u
CHARACTER(len=32)::fmt
WRITE(fmt,*)'(F5.2,',nt,'F5.2)'
!initial_condition
DO i = 1,nx
IF(i*dx>=0.5 .and. i*dx<=1) THEN
u(1,i)=2
ELSE
u(1,i)=1
ENDIF
ENDDO
!Finite_Difference
DO it=2,nt
DO i=2,nx-1
u(it,i)=u(it-1,i)-c*dt/dx*(u(it-1,i)-u(it-1,i-1))
ENDDO
ENDDO
!Write_File
OPEN(UNIT=200,FILE='tab2.txt',STATUS='REPLACE',ACTION='WRITE')
DO i=1,nx
WRITE(200,fmt)dble(i)*dx,u(:,i)
ENDDO
CLOSE(UNIT=200)
END PROGRAM odlc |
|
|
Back to top |
|
 |
|
|
You cannot post new topics in this forum You cannot reply to topics in this forum You cannot edit your posts in this forum You cannot delete your posts in this forum You cannot vote in polls in this forum
|
Powered by phpBB © 2001, 2005 phpBB Group
|