From 3368f302df4e43d5d49e10f2bca665f07d94fed1 Mon Sep 17 00:00:00 2001 From: bangerth Date: Thu, 2 Dec 2010 16:42:00 +0000 Subject: [PATCH] Use better gravity model as described in the introduction. Allow graphical output also for more processors but more infrequently. git-svn-id: https://svn.dealii.org/trunk@22904 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/examples/step-32/step-32.cc | 18 +++++++++--------- 1 file changed, 9 insertions(+), 9 deletions(-) diff --git a/deal.II/examples/step-32/step-32.cc b/deal.II/examples/step-32/step-32.cc index 0330e8d547..776417ae48 100644 --- a/deal.II/examples/step-32/step-32.cc +++ b/deal.II/examples/step-32/step-32.cc @@ -109,6 +109,9 @@ namespace EquationData const double T0 = 4000+273; const double T1 = 700+273; + const double g0 = 10.7; + const double g1 = 9.81; + const double year_in_seconds = 60*60*24*365.2425; const double end_time = 1e9 * year_in_seconds; @@ -126,7 +129,8 @@ namespace EquationData template Tensor<1,dim> gravity_vector (const Point &p) { - return -9.81 * p / R1; + const double r = p.norm(); + return -(1.245e-6 * r + 7.714e13/r/r) * p / p.norm(); } @@ -156,15 +160,11 @@ namespace EquationData const double s = (r-R0)/h; // see http://www.wolframalpha.com/input/?i=plot+(sqrt(x^2%2By^2)*0.95%2B0.05*sin(6*atan2(x,y))),+x%3D-1+to+1,+y%3D-1+to+1 - double s_mod = s*0.95 + 0.05*sin(6.0*atan2(p(0),p(1))); + const double s_mod = s*0.95 + 0.05*sin(6.0*atan2(p(0),p(1))); //alternative: http://www.wolframalpha.com/input/?i=plot+atan((sqrt(x^2%2By^2)*0.95%2B0.05*sin(6*atan2(x,y))-0.5)*10)/pi%2B0.5,+x%3D-1+to+1,+y%3D-1+to+1 // s_mod = atan((s_mod-0.5)*10.0)/dealii::numbers::PI+0.5; - return T1+(T0-T1)*(1.0-s_mod); - - //old: -// return T1+(T0-T1)*((1-s)*(1-s)); -// return T1+(T0-T1)*((1-s)); + return T0*(1.0-s_mod) + T1*s_mod; } @@ -3392,8 +3392,8 @@ void BoussinesqFlowProblem::run (unsigned int ref) if ((timestep_number > 0) && (timestep_number % 10 == 0)) refine_mesh (initial_refinement + n_pre_refinement_steps); - if (timestep_number % 10 == 0 && - Utilities::System::get_n_mpi_processes(MPI_COMM_WORLD) <= 10) + if (timestep_number % 50 == 0 && + Utilities::System::get_n_mpi_processes(MPI_COMM_WORLD) <= 100) output_results (); time += time_step; -- 2.39.5