# problem in writing a solver for a 3-step runge kutta(incompressible flow)

 User Name Remember Me Password
 Register Blogs Members List Search Today's Posts Mark Forums Read

 LinkBack Thread Tools Display Modes
 February 24, 2014, 14:03 problem in writing a solver for a 3-step runge kutta(incompressible flow) #1 Member   Join Date: Feb 2012 Posts: 49 Rep Power: 6 I'm implementing a 3-step runge kutta descritization for solving incompressible flow. The results are completely wrong. I would appreciate if you have a look at the code below and tell me if you see any problem : Code: ```while (runTime.loop()) { double alpha,gamma,zeta; volVectorField Unew(U), UOld2(U), UTmp(U) ; surfaceScalarField phiOld2(phi), phiTmp(phi); volScalarField pPhi(p); for (int i = 1 ; i<=3 ; ++i) { switch (i) { case 1 : {UOld2=U ; phiOld2=phi ; alpha=4./15. ; gamma= 8./15. ; zeta =0.; break ; } case 2 : {UTmp=U ; phiTmp=phi ; alpha=1./15. ; gamma= 5./12. ; zeta =-17./60.; break ; } case 3 : {UOld2=UTmp ; phiOld2=phiTmp ; alpha=1./6. ; gamma= 3./4. ; zeta =-5./12.; break ; } default : break ; } adjustPhi(phi, U, p); U = U + runTime.deltaT() * ( + 2*alpha*fvc::laplacian(nu,U) - 2*alpha*fvc::grad(p) - gamma*fvc::div(phi, U) - zeta*fvc::div(phiOld2, UOld2) ); solve(alpha*runTime.deltaT()*fvm::laplacian(nu,Unew) - fvm::Sp(1.,Unew) == //(alpha*runTime.deltaT()) == (-1.)*(U) + alpha*runTime.deltaT()*fvc::laplacian(nu,U) ); solve (fvm::laplacian(pPhi) == fvc::div(U)/(2.*alpha*runTime.deltaT())); U = Unew - (2.*alpha*runTime.deltaT()*fvc::grad(pPhi)); p += pPhi - alpha*runTime.deltaT()*nu*(fvc::laplacian(pPhi)); adjustPhi(phi, U, p); U.correctBoundaryConditions(); } runTime.write(); } }```

 March 4, 2014, 12:19 #2 Member   Join Date: Feb 2012 Posts: 49 Rep Power: 6 I make my code shorter, maybe someone has any idea regarding any problem on this kind of implementation : Code: ```while (runTime.loop()) { for (int i = 1 ; i<=3 ; ++i) { U = U + runTime.deltaT() * ( + 2*alpha*fvc::laplacian(nu,U) - 2*alpha*fvc::grad(p) - gamma*fvc::div(phi, U) - zeta*fvc::div(phiOld2, UOld2) ); solve(alpha*runTime.deltaT()*fvm::laplacian(nu,Unew) - fvm::Sp(1.,Unew) == //(alpha*runTime.deltaT()) == (-1.)*(U) + alpha*runTime.deltaT()*fvc::laplacian(nu,U) ); solve (fvm::laplacian(pPhi) == fvc::div(U)/(2.*alpha*runTime.deltaT()));// pPhi is a pseudo pressure without physical meaning U = Unew - (2.*alpha*runTime.deltaT()*fvc::grad(pPhi)); p += pPhi - alpha*runTime.deltaT()*nu*(fvc::laplacian(pPhi)); adjustPhi(phi, U, p); U.correctBoundaryConditions(); } }```

 Thread Tools Display Modes Linear Mode

 Posting Rules You may not post new threads You may not post replies You may not post attachments You may not edit your posts BB code is On Smilies are On [IMG] code is On HTML code is OffTrackbacks are On Pingbacks are On Refbacks are On Forum Rules

 Similar Threads Thread Thread Starter Forum Replies Last Post Lord Kelvin OpenFOAM Running, Solving & CFD 8 March 28, 2016 11:08 xiuying OpenFOAM Running, Solving & CFD 8 August 27, 2013 15:33 rae CFX 3 June 6, 2013 06:49 roo FLUENT 3 June 4, 2013 14:53 lzgwhy FLUENT 0 August 19, 2009 22:44

All times are GMT -4. The time now is 17:00.

 Contact Us - CFD Online - Top