Hi guys

I am making an compressible euler code using FVM on structured grids (i am going unstructured in the future,but i take things step by step.).I have made an explicit version of the code,and now i am working on the implicit version.

I am using first order Roe scheme for temporal discretisation and 2nd order MUSCL Roe/Van leer scheme for spatial.I solve the system using DADI.The problem is that the implicit version doesn't converge even with cf lower than one(the explicit code runs at larger cfls!!!)

It's obvious that i am doing sth wrong...So basically i am asking for any ideas of what is going wrong..

I have tested the implicit code in 1-d shock tube and it converges with max CFL 5.

Is there a chance that additional dissipation is needed even though i use upwind schemes??

Can it be that adi si causing the instabilities??

Thank you very much.

PS. I know that i am giving little information for you to point me the errors.I am just asking for any ideas of what might be wrong.People with greater experience might have in mind something

I am putting some code here,if anyone hcan help

!!This is the I-direction sweep(the same is for J)

Subroutine sweepi

use Variables

use Param

use Matr

use grid

Implicit None

real,allocatable :: VARP(:,:,:,,VARM(:,:,:,

real :: UV(3),E,ATMP(4,4),NLOC(3),UB(3),TMP(4,4),TMP1(4,4) ,ATMP1(4,4)

integer :: i,j

real :: P1,P2,c1,c2,dens1,dens2,Vt,U1(3),U2(3),p,c,dens

allocate (VARP( 1:NMAXI,1:NMAXI,NMAXK,5),VARM(1:NMAXI,1:NMAXJ,NMAX K,5))

UB=0

Vt=0

NLOC=0 ! First order temporal

do i=1,NMAXI-1

do j=1,NMAXJ-1

VARP(I,J,1,1:5) = UP(I ,J,1,1:5)

VARM(I,J,1,1:5) = UP(I-1 ,J,1,1:5)

NLOC(1:2) = NORMALI(I ,J,1:2)

!JACIP

dens = VARP(I,J,1,1)

UV(1:3) = VARP(I,J,1,2:4) / VARP(I,J,1,1)

E = VARP(I,J,1,5 ) / VARP(I,J,1,1)

P = (GAM-1)*(VARP(I,J,1,5) -0.5*VARP(I,J,1,1)*(UV(1)**2+UV(2)**2))

c = sqrt(GAM*p/VARP(I,J,1,1))

call conv_jacob2(nloc,dens,UV,E,ATMP)

JACIP(i,j,1:4,1:4) = AREA(I,J,1)*ATMP (1:4,1:4)

!

!JACIM

dens = VARM(I,J,1,1)

UV(1:3) = VARM(I,J,1,2:4) / VARM(I,J,1,1)

E = VARM(I,J,1,5 ) / VARM(I,J,1,1)

P = (GAM-1)*(VARM(I,J,1,5) -0.5*VARM(I,J,1,1)*(UV(1)**2+UV(2)**2))

c = sqrt(GAM*p/VARM(I,J,1,1))

call conv_jacob2(nloc,dens,UV,E,ATMP)

JACIM(i,j,1:4,1:4) = AREA(I,J,1)*ATMP (1:4,1:4)

!

!ROE

dens1 =VARP(I,J,1,1)

dens2 =VARM(I,J,1,1)

U1(1:2) =VARP(I,J,1,2:3)/dens1

U2(1:2) =VARM(I,J,1,2:3)/dens2

P1 = (GAM-1)*(VARP(i,j,1,5) -0.5*dens1*(U1(1)**2+U1(2)**2))

P2 = (GAM-1)*(VARM(i,j,1,5) -0.5*dens2*(U2(1)**2+U2(2)**2))

c1 = sqrt(GAM*p1/dens1)

c2 = sqrt(GAM*p2/dens2)

call roejac2 (dens1,dens2,c1,c2 ,U1,U2,P1,P2,NLOC,Vt,ATMP)

JACIROE(i,j,1:4,1:4) =AREA(i,j,1)*ATMP(1:4,1:4)

enddo

enddo

deallocate (VARP,VARM)

End Subroutine sweepi

!the construction of the implicit operator

Subroutine impl2(i,j,nodi,nodj)

use Variables

use Param

use Matr

use grid

Implicit None

integer,intent(in) :: i,j,nodi,nodj

real :: AI(4,4),II(1:4,1:4),W(1:4,1:4)

AI=0

II=0

II(1,1)=1

II(2,2)=1

II(3,3)=1

II(4,4)=1

AI(1,1)=VOLUME(i,j)/DTT(i,j)

AI(2,2)=VOLUME(i,j)/DTT(i,j)

AI(3,3)=VOLUME(i,j)/DTT(i,j)

AI(4,4)=VOLUME(i,j)/DTT(i,j)

LOWERI(NODI,1,1:4,1:4 ) = 0.5*( JACIM(i,j,1:4,1:4) - JACIROE(i,j,1:4,1:4) )

LOWERJ(NODJ,1,1:4,1:4 ) = 0.5*( JACJM(i,j,1:4,1:4) - JACJROE(i,j,1:4,1:4) )

UPPERI(NODI,1,1:4,1:4 ) =0.5*(-JACIP(i+1,j ,1:4,1:4) - JACIROE(i+1,j ,1:4,1:4))

UPPERJ(NODJ,1,1:4,1:4 ) =0.5*(-JACJP(i ,j+1,1:4,1:4) - JACJROE(i ,j+1,1:4,1:4))

DIAGI(NODI,1,1:4,1:4)=0.5*(-JACIM(i+1,j,1:4,1:4)+JACIP(i,j,1:4,1:4)+JACIROE(i ,j,1:4,1:4) +JACIROE(i+1,j,1:4,1:4) ) +AI(1:4,1:4)

DIAGJ(NODJ,1,1:4,1:4)=0.5*(-JACJM(i,j+1,1:4,1:4)+JACJP(i,j,1:4,1:4)+JACJROE(i, j,1:4,1:4)+JACJROE(i,j+1,1:4,1:4))+AI(1:4,1:4)

end Subroutine

[Nishu]

Hi,

Can you please send me the 1D code for shock tube problem at babuu.nishu@gmail.com? I need it for my research.

Thanks alot
Nishant