Hello David, Thank you for the answer. I know Step-35 doesn't implement velocity correction, which
is still confusing me. Without correcting the velocity, it is not divergence free, does it make any sense? Best, Jiaqi 2016-08-02 16:37 GMT-04:00 David Wells <drwe...@vt.edu>: > Hi Jiaqi, > > I have to admit that it has been about a year since I last really looked > at step-35 and I do not know the answer off the top of my head for this. I > saw that Daniel responded to this question elsewhere; was that enough to > figure things out? > > Thanks, > David Wells > > > On Friday, July 29, 2016 at 9:19:29 PM UTC-4, Jiaqi ZHANG wrote: >> >> Hello David Wells, >> >> Sorry to bother you, but I have a problem about Step-35. I have posted >> the question >> >> a few days ago, and no one answered me, so I have to search in the >> mailing list and found >> >> that you have been using it a lot. I was wondering if you could help me, >> and the following is >> >> my question: >> >> I read the references listed there and some other papers about projection >> methods, >> >> and know that after obtaining the intermediate values for velocity, say >> v_n, we need to >> >> correct the velocity by u_n = v_n - grad(2dt/3 phi_n), which is mentioned >> in the introduction of Step-35, >> >> but it seems that Step-35 only corrects the pressure, and I cannot find >> anything in the code that implements the above equation. >> >> However, when I correct the velocity at every time step, the algorithm >> turns out to be divergent (iterating more than 1000 times at some time >> step). >> >> I am very confused, did I misunderstand the algorithm or the code? >> >> Thanks in advance. >> >> >> Best, >> Jiaqi >> >> >> >> 在 2015年3月17日星期二 UTC-4下午8:14:21，David Wells写道： >>> >>> Hi Joaquin, >>> >>> I used step-35 a lot, so I may be able to assist you here. >>> >>> For the velocity inflow, are you sure that what you have will be zero >>> flow at the edges? The XY slices don't look quite circular to me, so this >>> may be a problem if you refine the mesh. Since you are using GMSH, it may >>> be easiest to use GMSH's extrusion abilities to extrude a circle along the >>> z axis to solve this. >>> >>> I believe you handled the pressure correctly (as long as the domain >>> starts at z = 0 and ends at z = 10). >>> >>> At least as of a few months ago, step-35 could not run in 3D due to a >>> lack of a compressed sparsity pattern. I fixed this in my own copy of >>> step-35 (see http://www.github.com/drwells/step-35a) (yeah, I know I >>> should submit a patch...). I tweaked the preconditioners to improve >>> performance, so my version might be useful to look at. >>> >>> Another thing you might want to look at is how Abner implemented labels >>> for the different boundaries. According to my experiments (see >>> https://github.com/drwells/dealii-step-35a/blob/master/cylinderextruded.geo >>> ): >>> 1 is the no-slip edges >>> 2 is the inflow >>> 3 is the outflow >>> 4 is the cylinder (may not apply for you) >>> >>> If your mesh file does not label things correctly then things will get >>> weird. I hope that this is helpful. >>> >>> Thanks, >>> David Wells >>> >>> On Tuesday, March 17, 2015 at 6:08:02 PM UTC-4, Joaquin wrote: >>>> >>>> Hi, >>>> >>>> I am trying to run step-35 for the unsteady incompressible viscous flow >>>> in a pipe, for 3D case. The mesh file I used is attached. The graph of >>>> results >>>> is attached. The expected results for any time should be paraboloid >>>> (or parabolic in 2D), but, I can not get that. The changes I did are: >>>> >>>> >>>> >>>> *1. namespace EquationData: * template <int dim> >>>> double Velocity<dim>::value (const Point<dim> &p, >>>> const unsigned int) const >>>> { >>>> if (this->comp == 0) >>>> { >>>> const double Um = 1.5; >>>> const double R = 2.; >>>> return Um*(1 - (p(0)*p(0)+p(1)*p(1))/(2*R*R)); >>>> } >>>> else >>>> return 0.; >>>> } >>>> >>>> *Through the pipe pressure drop **is linear (**I do not how to make >>>> changes here):* >>>> >>>> template <int dim> >>>> double Pressure<dim>::value (const Point<dim> &p, >>>> const unsigned int) const >>>> { >>>> return 10.-p(2); >>>> } >>>> >>>> >>>> *2. * >>>> template <int dim> >>>> void >>>> NavierStokesProjection<dim>:: >>>> create_triangulation_and_dofs (const unsigned int n_refines) >>>> { >>>> GridIn<dim> grid_in; >>>> grid_in.attach_triangulation (triangulation); >>>> >>>> { >>>> std::string filename = "*cilindro.msh*"; >>>> std::ifstream file (filename.c_str()); >>>> Assert (file, ExcFileNotOpen (filename.c_str())); >>>> grid_in.read_msh (file); >>>> } >>>> >>>> * 3. * >>>> int main() >>>> { >>>> try >>>> { >>>> using namespace dealii; >>>> using namespace Stepw3; >>>> >>>> RunTimeParameters::Data_Storage data; >>>> data.read_data ("parameter-file.prm"); >>>> >>>> deallog.depth_console (data.verbose ? 2 : 0); >>>> >>>> NavierStokesProjection<3> test (data); >>>> test.run (data.verbose, data.output_interval); >>>> } >>>> . >>>> . >>>> . >>>> return 0; >>>> } >>>> 4. >>>> template <int dim> >>>> void NavierStokesProjection<dim>::assemble_vorticity (const bool >>>> reinit_prec) >>>> { >>>> Assert (dim == 3, ExcNotImplemented()); >>>> if (reinit_prec) >>>> . >>>> . >>>> >>>> Into the parameter file: >>>> - Set n_of_refines = 1 >>>> - Set max_iterations = 5000 # maximal number of iterations that GMRES >>>> must make >>>> >>>> Please, can someone give me some guidelines to achieve my goal?, I >>>> have one week to present this easy step (for me it is difficult). I have >>>> studied the projection step method, but, the details I can't understand >>>> very well yet. (The boundary conditions are: at r = 0 then v = vmax; r = R >>>> then v = 0, Initial condition: t = 0 then v = 0). >>>> >>>> thanks in advance, >>>> Joaquin >>>> >>>> -- > The deal.II project is located at http://www.dealii.org/ > For mailing list/forum options, see > https://groups.google.com/d/forum/dealii?hl=en > --- > You received this message because you are subscribed to the Google Groups > "deal.II User Group" group. > To unsubscribe from this group and stop receiving emails from it, send an > email to dealii+unsubscr...@googlegroups.com. > For more options, visit https://groups.google.com/d/optout. > -- Jiaqi Zhang Graduate student Department of Mathematics Virginia Tech -- The deal.II project is located at http://www.dealii.org/ For mailing list/forum options, see https://groups.google.com/d/forum/dealii?hl=en --- You received this message because you are subscribed to the Google Groups "deal.II User Group" group. To unsubscribe from this group and stop receiving emails from it, send an email to dealii+unsubscr...@googlegroups.com. For more options, visit https://groups.google.com/d/optout.
