From 8141de7d5de4ec7b6089eba4c1d0f54daa7b6aad Mon Sep 17 00:00:00 2001 From: Ian Hinder Date: Thu, 13 Oct 2011 02:39:35 +0200 Subject: Add Laplace equation example This uses Jacobi iteration to solve Laplace's equation in 2D on a square grid. --- Examples/Laplace.m | 68 +++++++++++++++++++++++++++++++++++++++ Examples/laplace.par | 91 ++++++++++++++++++++++++++++++++++++++++++++++++++++ 2 files changed, 159 insertions(+) create mode 100644 Examples/Laplace.m create mode 100644 Examples/laplace.par (limited to 'Examples') diff --git a/Examples/Laplace.m b/Examples/Laplace.m new file mode 100644 index 0000000..ad8a28f --- /dev/null +++ b/Examples/Laplace.m @@ -0,0 +1,68 @@ + +groups = {{"phi_group", {phi}}}; + +derivatives = +{ + PDstandard[i_] -> + StandardCenteredDifferenceOperator[1,fdOrder/2,i], + PDstandard[i_, i_] -> + StandardCenteredDifferenceOperator[2,fdOrder/2,i], + PDstandard[i_, j_] -> + StandardCenteredDifferenceOperator[1,fdOrder/2,i] StandardCenteredDifferenceOperator[1,fdOrder/2,j] +}; + +PD = PDstandard; + +initialCalc = +{ + Name -> "Laplace_initial", + Schedule -> {"AT INITIAL"}, + Where -> Interior, + Equations -> + { + phi -> phi0 Sum[4/(Pi n) Sin[n Pi x/Lx] Sinh[n Pi y/Lx]/Sinh[n Pi Ly/Lx], {n, 1, 1, 2}] + } +}; + +initialBoundaryCalc = +{ + Name -> "Laplace_initial_boundary", + Schedule -> {"AT INITIAL after Laplace_initial"}, + Where -> Boundary, + Equations -> + { + phi -> IfThen[Abs[y-Ly]<10^-10, phi0, 0] + } +}; + +evolveCalc = +{ + Name -> "Laplace_relax", + Schedule -> {"in MoL_CalcRHS", "AT ANALYSIS"}, + Equations -> + { + dot[phi] -> mu Euc[ui,uj] PD[phi,li,lj] + } +}; + +boundaryCalc = +{ + Name -> "Laplace_boundary", + Schedule -> {"in MoL_RHSBoundaries", "AT ANALYSIS"}, + Where -> Boundary, + Equations -> + { + dot[phi] -> 0 + } +}; + + +CreateKrancThornTT[ + groups, ".", + "Laplace", + Calculations -> {initialCalc, initialBoundaryCalc, evolveCalc, boundaryCalc}, + PartialDerivatives -> derivatives, + ZeroDimensions -> {3}, + RealParameters -> {Lx,Ly,phi0,mu}, + IntParameters -> {{Name -> fdOrder, Default -> 2, AllowedValues -> {2, 4}}}, + DeclaredGroups -> {"phi_group"}]; diff --git a/Examples/laplace.par b/Examples/laplace.par new file mode 100644 index 0000000..a63db80 --- /dev/null +++ b/Examples/laplace.par @@ -0,0 +1,91 @@ + +Cactus::cctk_final_time = 1 +Cactus::terminate = "time" + +ActiveThorns = "IOUtil Carpet CarpetLib CarpetSlab CoordBase CoordBase SymBase CartGrid3D Slab CarpetIOBasic CarpetIOASCII CarpetIOScalar Laplace Time MoL Boundary GenericFD CarpetReduce LoopControl CarpetIOHDF5 HDF5" + +CoordBase::domainsize = minmax + +CoordBase::boundary_size_x_lower = 1 +CoordBase::boundary_size_y_lower = 1 +CoordBase::boundary_size_z_lower = 0 +CoordBase::boundary_shiftout_x_lower = 0 +CoordBase::boundary_shiftout_y_lower = 0 +CoordBase::boundary_shiftout_z_lower = 1 + +CoordBase::boundary_size_x_upper = 1 +CoordBase::boundary_size_y_upper = 1 +CoordBase::boundary_size_z_upper = 0 +CoordBase::boundary_shiftout_x_upper = 0 +CoordBase::boundary_shiftout_y_upper = 0 +CoordBase::boundary_shiftout_z_upper = 1 + +CartGrid3D::type = "coordbase" +CartGrid3D::domain = "full" +CartGrid3D::avoid_origin = "no" + +CoordBase::xmin = 0 +CoordBase::ymin = 0 +CoordBase::zmin = 0 +CoordBase::xmax = 1 +CoordBase::ymax = 1 +CoordBase::zmax = 0 +CoordBase::dx = 0.02 +CoordBase::dy = 0.02 +CoordBase::dz = 0.02 + +driver::ghost_size_x = 1 +driver::ghost_size_y = 1 +driver::ghost_size_z = 0 + +Carpet::domain_from_coordbase = "yes" +Carpet::poison_new_timelevels = "yes" +Carpet::check_for_poison = "no" + +Carpet::max_refinement_levels = 1 +Carpet::prolongation_order_space = 3 +Carpet::prolongation_order_time = 2 +Carpet::init_each_timelevel = no +Carpet::print_timestats_every = 0 + +Time::timestep = 0.0001 +Time::timestep_method = "given" + +MethodOfLines::ode_method = "generic" +MethodOfLines::generic_type = "RK" +MethodOfLines::MoL_Intermediate_Steps = 1 +MethodOfLines::MoL_Num_Scratch_Levels = 1 + +# MethodOfLines::ode_method = "RK2" +# MethodOfLines::MoL_NaN_Check = "no" +# MethodOfLines::initial_data_is_crap = "yes" +# MethodOfLines::MoL_Intermediate_Steps = 2 +# MethodOfLines::MoL_Num_Scratch_Levels = 1 + +IO::out_dir = $parfile + +IOBasic::outInfo_every = 100 +IOBasic::outInfo_reductions = "minimum maximum" +IOBasic::outInfo_vars = " + Carpet::physical_time_per_hour + Laplace::phi + Laplace::phirhs +" + +IOScalar::outScalar_every = 10 +IOScalar::outScalar_reductions = "norm2" +IOScalar::outScalar_vars = "Laplace::phirhs" + +Laplace::Lx = 1 +Laplace::Ly = 1 +Laplace::phi0 = 1 +Laplace::mu = 1 + +IOHDF5::out_every = 100 +IOHDF5::out_vars = " + Laplace::phi + Laplace::phirhs +" + +# TimerReport::out_every = 0 +# TimerReport::n_top_timers = 50 -- cgit v1.2.3