From 53a97e2f7fc607043b45f333a876bfeece7a004b Mon Sep 17 00:00:00 2001 From: Lynn Munday Date: Thu, 14 Sep 2023 18:16:56 -0600 Subject: [PATCH] Update python script to solve Gringarten solution for finite spaced fractures and update input files for this case. I also removed examples that used AMR because they were not useful. closes #137 --- examples/gringarten_multiapp/frect1.i | 120 +++--- examples/gringarten_multiapp/frect1_amr.i | 386 ------------------ .../gringartenbrickdomain100by100by10.e | Bin 107184 -> 0 bytes examples/gringarten_multiapp/mrect1.i | 98 +++-- examples/gringarten_multiapp/mrect1_amr.i | 265 ------------ examples/gringarten_multiapp/plot_compare.py | 142 +++++-- .../production_single_fracture_amr.xyz | 1 - .../gringarten_multiapp/rect_100_10_medium.e | Bin 14484 -> 0 bytes examples/gringarten_multiapp/tests | 7 +- 9 files changed, 250 insertions(+), 769 deletions(-) delete mode 100644 examples/gringarten_multiapp/frect1_amr.i delete mode 100644 examples/gringarten_multiapp/gringartenbrickdomain100by100by10.e delete mode 100644 examples/gringarten_multiapp/mrect1_amr.i delete mode 100644 examples/gringarten_multiapp/production_single_fracture_amr.xyz delete mode 100644 examples/gringarten_multiapp/rect_100_10_medium.e diff --git a/examples/gringarten_multiapp/frect1.i b/examples/gringarten_multiapp/frect1.i index 8e8680dd..855c0ceb 100644 --- a/examples/gringarten_multiapp/frect1.i +++ b/examples/gringarten_multiapp/frect1.i @@ -1,17 +1,29 @@ # Cold water injection into one side of the fracture network, and production from the other side -injection_rate = 0.1 #25 # kg/s -endTime = 3.16e8 +injection_rate = 0.1 # kg/s +endTime = 3e8 [Mesh] uniform_refine = 0 - [single_frac] - type = FileMeshGenerator - file = 'rect_100_10_medium.e' + [generate] + type = GeneratedMeshGenerator + dim = 2 + nx = 10 + xmin = -5 + xmax = 5 + ny = 50 + ymin = -50 + ymax = 50 + [] + [rotate] + type = TransformGenerator + input = generate + transform = ROTATE + vector_value = '90 90 90' [] [injection_node] type = BoundingBoxNodeSetGenerator - input = single_frac - bottom_left = '-5 -55 -5' - top_right = '5 -45 5' + input = rotate + bottom_left = '-0.05 -51 -5.1' + top_right = '0.05 -49.9 5.1' new_boundary = injection_node [] [] @@ -45,8 +57,9 @@ endTime = 3.16e8 coupling_type = ThermoHydro porepressure = frac_P temperature = frac_T - fp = water + fp = the_simple_fluid #water pressure_unit = Pa + # stabilization = full [] [Kernels] @@ -129,10 +142,10 @@ endTime = 3.16e8 [heat_transfer_coefficient_auxk] type = ParsedAux variable = heat_transfer_coefficient - args = 'enclosing_element_normal_length enclosing_element_normal_thermal_cond' + coupled_variables = 'enclosing_element_normal_length enclosing_element_normal_thermal_cond' constant_names = h_s constant_expressions = 1E3 #This is the value being assigned to h_s. Should be much bigger than thermal_conductivity / L ~ 1 - function = 'if(enclosing_element_normal_length = 0, 0, h_s * enclosing_element_normal_thermal_cond * 2 * enclosing_element_normal_length / (h_s * enclosing_element_normal_length * enclosing_element_normal_length + enclosing_element_normal_thermal_cond * 2 * enclosing_element_normal_length))' + expression = 'if(enclosing_element_normal_length = 0, 0, h_s * enclosing_element_normal_thermal_cond * 2 * enclosing_element_normal_length / (h_s * enclosing_element_normal_length * enclosing_element_normal_length + enclosing_element_normal_thermal_cond * 2 * enclosing_element_normal_length))' [] [insitu_pp] type = FunctionAux @@ -219,11 +232,17 @@ endTime = 3.16e8 [] [FluidProperties] + [the_simple_fluid] + type = SimpleFluidProperties + bulk_modulus = 2E9 + viscosity = 1.0E-3 + density0 = 1000.0 + [] [true_water] type = Water97FluidProperties [] [water] - type = TabulatedFluidProperties + type = TabulatedBicubicFluidProperties fp = true_water temperature_min = 275 # K temperature_max = 600 @@ -233,27 +252,21 @@ endTime = 3.16e8 [] [Materials] - [porosity] - type = PorousFlowPorosityLinear - porosity_ref = 1E-4 - P_ref = insitu_pp - P_coeff = 3e-10 - porosity_min = 1E-5 - [] - [permeability] - type = PorousFlowPermeabilityKozenyCarman - k0 = 1E-15 - poroperm_function = kozeny_carman_phi0 - m = 0 - n = 3 - phi0 = 1E-4 - [] - [internal_energy] + [porosity_frac] + type = PorousFlowPorosityConst + porosity = 0.01 + [] + [permeability_frac] + type = PorousFlowPermeabilityConst + permeability = '1e-12 0 0 0 1e-12 0 0 0 1e-12' + [] + + [internal_energy_frac] type = PorousFlowMatrixInternalEnergy density = 2700 specific_heat_capacity = 0 [] - [aq_thermal_conductivity] + [aq_thermal_conductivity_frac] type = PorousFlowThermalConductivityIdeal dry_thermal_conductivity = '0.6E-4 0 0 0 0.6E-4 0 0 0 0.6E-4' [] @@ -262,17 +275,17 @@ endTime = 3.16e8 [Functions] [kg_rate] type = ParsedFunction - vals = 'dt kg_out' - vars = 'dt kg_out' - value = 'kg_out/dt' + symbol_values = 'dt kg_out' + symbol_names = 'dt kg_out' + expression = 'kg_out/dt' [] [insitu_pp] type = ParsedFunction - value = '9.81*1000*(3000 - z)' + expression = '9.81*1000*(3000 - z)' [] [insitu_T] type = ParsedFunction - value = '363' + expression = '363' [] [] @@ -284,19 +297,23 @@ endTime = 3.16e8 [kg_out] type = PorousFlowPlotQuantity uo = kg_out_uo + outputs = none [] [kg_per_s] type = FunctionValuePostprocessor function = kg_rate + outputs = none [] [J_out] type = PorousFlowPlotQuantity uo = J_out_uo + outputs = none [] [TK_in] type = PointValue variable = frac_T point = '0 -50 0' + outputs = none [] [TK_out] type = PointValue @@ -307,36 +324,46 @@ endTime = 3.16e8 type = PointValue variable = frac_P point = '0 50 0' + outputs = none [] [P_in] type = PointValue variable = frac_P point = '0 -50 0' + outputs = none [] [] [VectorPostprocessors] [heat_transfer_rate] type = NodalValueSampler - outputs = none sort_by = id variable = joules_per_s + outputs = none [] [] [Preconditioning] - [./superlu] + [superlu] type = SMP full = true petsc_options_iname = '-ksp_type -pc_type -pc_factor_mat_solver_package' petsc_options_value = 'gmres lu superlu_dist' - [../] + [] [] +[Functions] + [dts] + type = ParsedFunction + expression = if(t<1e6,t*t,1e6) + [] +[] [Executioner] type = Transient solve_type = NEWTON - dt = 0.1e7 + # dt = 1e6 + dtmin = 1 + dtmax = 1e6 end_time = ${endTime} line_search = 'none' automatic_scaling = true @@ -346,21 +373,14 @@ endTime = 3.16e8 nl_max_its = 20 nl_rel_tol = 5e-04 nl_abs_tol = 1e-09 + [TimeStepper] + type = FunctionDT + function = dts + min_dt = 2 + [] [] [Outputs] print_linear_residuals = false - exodus = false csv = true - # [fracCSV] - # type = CSV - # sync_times = '100 200 300 400 500 600 700 800 900 - # 1000 2000 3000 4000 5000 6000 7000 8000 9000 - # 1000e1 2000e1 3000e1 4000e1 5000e1 6000e1 7000e1 8000e1 9000e1 - # 1000e2 2000e2 3000e2 4000e2 5000e2 6000e2 7000e2 8000e2 9000e2 - # 1000e3 2000e3 3000e3 4000e3 5000e3 6000e3 7000e3 8000e3 9000e3 - # 1000e4 2000e4 3000e4 4000e4 5000e4 6000e4 7000e4 8000e4 9000e4 - # 1000e5 2000e5 3000e5 4000e5 5000e5 6000e5 7000e5 8000e5 9000e5' - # sync_only = true - # [] [] diff --git a/examples/gringarten_multiapp/frect1_amr.i b/examples/gringarten_multiapp/frect1_amr.i deleted file mode 100644 index b5ea5cae..00000000 --- a/examples/gringarten_multiapp/frect1_amr.i +++ /dev/null @@ -1,386 +0,0 @@ -# Cold water injection into one side of the fracture network, and production from the other side -injection_rate = 0.1 #25 # kg/s -endTime = 3.16e8 - - -[Mesh] - uniform_refine = 0 - [generate] - type = GeneratedMeshGenerator - dim = 2 - nx = 10 - xmin = -5 - xmax = 5 - ny = 50 - ymin = -50 - ymax = 50 - [] - [./rotate] - type = TransformGenerator - input = generate - transform = ROTATE - vector_value = '90 90 90' - [] - [./offset] - type = TransformGenerator - input = rotate - transform = TRANSLATE - vector_value = '0.1 0 0' - [] - [injection_node] - type = BoundingBoxNodeSetGenerator - input = offset - bottom_left = '0.05 -51 -5.1' - top_right = '0.15 -49.9 5.1' - new_boundary = injection_node - [] -[] - -[GlobalParams] - PorousFlowDictator = dictator - gravity = '0 0 -9.81' -[] - -[Variables] - [frac_P] - [] - [frac_T] - [] -[] - -[ICs] - [frac_P] - type = FunctionIC - variable = frac_P - function = insitu_pp - [] - [frac_T] - type = FunctionIC - variable = frac_T - function = insitu_T - [] -[] - -[PorousFlowFullySaturated] - coupling_type = ThermoHydro - porepressure = frac_P - temperature = frac_T - fp = water - pressure_unit = Pa -[] - -[Kernels] - [toMatrix] - type = PorousFlowHeatMassTransfer - variable = frac_T - v = transferred_matrix_T - transfer_coefficient = heat_transfer_coefficient - save_in = joules_per_s - [] -[] - -[AuxVariables] - [transferred_matrix_T] - initial_condition = 363 - [] - [heat_transfer_coefficient] - family = MONOMIAL - order = CONSTANT - initial_condition = 0.0 - [] - [joules_per_s] - [] - [aperture] - family = MONOMIAL - order = CONSTANT - [] - [perm_times_app] - family = MONOMIAL - order = CONSTANT - [] - [density] - family = MONOMIAL - order = CONSTANT - [] - [viscosity] - family = MONOMIAL - order = CONSTANT - [] - [insitu_pp] - [] - [normal_dirn_x] - family = MONOMIAL - order = CONSTANT - [] - [normal_dirn_y] - family = MONOMIAL - order = CONSTANT - [] - [normal_dirn_z] - family = MONOMIAL - order = CONSTANT - [] - [enclosing_element_normal_length] - family = MONOMIAL - order = CONSTANT - [] - [enclosing_element_normal_thermal_cond] - family = MONOMIAL - order = CONSTANT - [] -[] - -[AuxKernels] - [normal_dirn_x_auxk] - type = PorousFlowElementNormal - variable = normal_dirn_x - component = x - [] - [normal_dirn_y] - type = PorousFlowElementNormal - variable = normal_dirn_y - component = y - [] - [normal_dirn_z] - type = PorousFlowElementNormal - variable = normal_dirn_z - component = z - [] - [heat_transfer_coefficient_auxk] - type = ParsedAux - variable = heat_transfer_coefficient - args = 'enclosing_element_normal_length enclosing_element_normal_thermal_cond' - constant_names = h_s - constant_expressions = 1E3 #This is the value being assigned to h_s. Should be much bigger than thermal_conductivity / L ~ 1 - function = 'if(enclosing_element_normal_length = 0, 0, h_s * enclosing_element_normal_thermal_cond * 2 * enclosing_element_normal_length / (h_s * enclosing_element_normal_length * enclosing_element_normal_length + enclosing_element_normal_thermal_cond * 2 * enclosing_element_normal_length))' - [] - [insitu_pp] - type = FunctionAux - execute_on = initial - variable = insitu_pp - function = insitu_pp - [] - [aperture] - type = PorousFlowPropertyAux - variable = aperture - property = porosity - [] - [perm_times_app] - type = PorousFlowPropertyAux - variable = perm_times_app - property = permeability - row = 0 - column = 0 - [] - [density] - type = PorousFlowPropertyAux - variable = density - property = density - phase = 0 - [] - [viscosity] - type = PorousFlowPropertyAux - variable = viscosity - property = viscosity - phase = 0 - [] -[] - -[BCs] - [inject_heat] - type = DirichletBC - boundary = injection_node - variable = frac_T - value = 303 - [] -[] - -[DiracKernels] - [inject_fluid] - type = PorousFlowPointSourceFromPostprocessor - mass_flux = ${injection_rate} - point = '0.1 -50 0' - variable = frac_P - [] - [withdraw_fluid] - type = PorousFlowPeacemanBorehole - SumQuantityUO = kg_out_uo - bottom_p_or_t = 30.6e6 - character = 1 - line_length = 1 - point_file = production_single_fracture_amr.xyz - unit_weight = '0 0 0' - fluid_phase = 0 - use_mobility = true - variable = frac_P - [] - [withdraw_heat] - type = PorousFlowPeacemanBorehole - SumQuantityUO = J_out_uo - bottom_p_or_t = 30.6e6 - character = 1 - line_length = 1 - point_file = production_single_fracture_amr.xyz - unit_weight = '0 0 0' - fluid_phase = 0 - use_mobility = true - use_enthalpy = true - variable = frac_T - [] -[] - -[UserObjects] - [kg_out_uo] - type = PorousFlowSumQuantity - [] - [J_out_uo] - type = PorousFlowSumQuantity - [] -[] - -[FluidProperties] - [true_water] - type = Water97FluidProperties - [] - [water] - type = TabulatedFluidProperties - fp = true_water - temperature_min = 275 # K - temperature_max = 600 - interpolated_properties = 'density viscosity enthalpy internal_energy' - fluid_property_file = water97_tabulated.csv - [] -[] - -[Materials] - [porosity] - type = PorousFlowPorosityLinear - porosity_ref = 1E-4 - P_ref = insitu_pp - P_coeff = 3e-10 - porosity_min = 1E-5 - [] - [permeability] - type = PorousFlowPermeabilityKozenyCarman - k0 = 1E-15 - poroperm_function = kozeny_carman_phi0 - m = 0 - n = 3 - phi0 = 1E-4 - [] - [internal_energy] - type = PorousFlowMatrixInternalEnergy - density = 2700 - specific_heat_capacity = 0 - [] - [aq_thermal_conductivity] - type = PorousFlowThermalConductivityIdeal - dry_thermal_conductivity = '0.6E-4 0 0 0 0.6E-4 0 0 0 0.6E-4' - [] -[] - -[Functions] - [kg_rate] - type = ParsedFunction - vals = 'dt kg_out' - vars = 'dt kg_out' - value = 'kg_out/dt' - [] - [insitu_pp] - type = ParsedFunction - value = '9.81*1000*(3000 - z)' - [] - [insitu_T] - type = ParsedFunction - value = '363' - [] -[] - -[Postprocessors] - [dt] - type = TimestepSize - outputs = 'none' - [] - [kg_out] - type = PorousFlowPlotQuantity - uo = kg_out_uo - [] - [kg_per_s] - type = FunctionValuePostprocessor - function = kg_rate - [] - [J_out] - type = PorousFlowPlotQuantity - uo = J_out_uo - [] - [TK_in] - type = PointValue - variable = frac_T - point = '0.1 -50 0' - [] - [TK_out] - type = PointValue - variable = frac_T - point = '0.1 50 0' - [] - [P_in] - type = PointValue - variable = frac_P - point = '0.1 -50 0' - [] - [P_out] - type = PointValue - variable = frac_P - point = '0.1 50 0' - [] -[] - -[VectorPostprocessors] - [heat_transfer_rate] - type = NodalValueSampler - outputs = none - sort_by = id - variable = joules_per_s - [] -[] - -[Preconditioning] - [./superlu] - type = SMP - full = true - petsc_options_iname = '-ksp_type -pc_type -pc_factor_mat_solver_package' - petsc_options_value = 'gmres lu superlu_dist' - [../] -[] - -[Executioner] - type = Transient - solve_type = NEWTON - dt = 0.1e7 - end_time = ${endTime} - line_search = 'none' - automatic_scaling = true - l_max_its = 20 - l_tol = 8e-3 - nl_forced_its = 1 - nl_max_its = 20 - nl_rel_tol = 5e-04 - nl_abs_tol = 1e-09 -[] - -[Outputs] - print_linear_residuals = false - exodus = false - csv = true - # [fracCSV] - # type = CSV - # sync_times = '100 200 300 400 500 600 700 800 900 - # 1000 2000 3000 4000 5000 6000 7000 8000 9000 - # 1000e1 2000e1 3000e1 4000e1 5000e1 6000e1 7000e1 8000e1 9000e1 - # 1000e2 2000e2 3000e2 4000e2 5000e2 6000e2 7000e2 8000e2 9000e2 - # 1000e3 2000e3 3000e3 4000e3 5000e3 6000e3 7000e3 8000e3 9000e3 - # 1000e4 2000e4 3000e4 4000e4 5000e4 6000e4 7000e4 8000e4 9000e4 - # 1000e5 2000e5 3000e5 4000e5 5000e5 6000e5 7000e5 8000e5 9000e5' - # sync_only = true - # [] -[] diff --git a/examples/gringarten_multiapp/gringartenbrickdomain100by100by10.e b/examples/gringarten_multiapp/gringartenbrickdomain100by100by10.e deleted file mode 100644 index 91ce521e2f22ed901a8eac18c1c49c124e95fce0..0000000000000000000000000000000000000000 GIT binary patch literal 0 HcmV?d00001 literal 107184 zcmeFZcbF8#_r5(dyDST$pcpWMISVK%BFZS3bIu|n0!k84F=Nh(IV%bxX3SZ}oD~%n z6%}(<%!si2-cQ##y}bqVS@idP|6;wa+tc;j)nR(Nrn+ZVyZ6|#kPH7brz`-~1Nsc= zIjGmbKDk_8rrGfR1N-zGHoVUewZcD52OT-E=OO(Es#O8a)E+eWkUql_U;QklT3Md~ zedtuV6z%4r-Sgl9hbNu1@r5kwKS*&$<+c=_E&!K&Kt8a3E zm4^;KawyA{&|JqpZ1~XrgZfJ6mWLM56c*?;q<_z&`V1Y`fAAohxh9a`EO+du)TU+f zHir%v+-o>n>p5g_|3Slh9zA&IAw7rnKbCWo>`$*j$JOl6{sTBiL1dhK>w&#S^dG46 z!#VBA0rnl>o{BJv!ZAH{(cg)1Gc|nD=p;mCymy>MQI!U%|B&p_oElN z{+@lumFm~A)*e5{diDF=T*lvG>zDJs{aebc-;nq1-&AJ(n{EI9b^GHRSl{r|MR;2a>v)u>+|iG+rM8=<+ksyuiX0n{{2_q(DwcK%5DGu%9nfn zzWs9JV_to^?a%A$A6w@7^8J?^Kfe8P`}gaoq3xG@{bL&V{ApZKlkOXe?Ncaw(sX}TxrwREhm5R$uCx?Tm3z1d;U`c zKJMpjuRqxOzc#Ra_r5=G`*z~v`>o9O{B`+uy#8$cf3)hk_jm8rkK6nEad<7~{W!~Q z*SGK6@qXU+8t1-VT)m&~*Vp&s_T%*2YwXKy*Uyvp_kNyx_2cmEdyU_x_&EK2_U*)N z-;P)J$NR_C_g}7``?}Y&!R$UueRa*aQmA928s~9)e3||?S?>Mgw(sYbX#H|+ zU*Go^w|#r=z3uhJXg_13wr|h-yZ7Vw^}QV*r}vN7_qJDG->c_0N9W6ZT)m&K7kkgW zHgsJ6I($3v>x|bgcbpcHDqs1AzSaX)YOw|?IC?ZoYY)~|2W|FEbXukZa2w|?$@J@>vp@8@klKHtBu z=daJbujg%F-`DfDd*6<)=k2)8Yk%)ouEyJq_v`D${k&bSzrU_>uhZA_?fZH2?fQD| zeS5y1xBY$c_56MGe(t@$dvE)C-uC{!|G0nLj*l<)-uBn$ZQozN((0ECoOa*wQ>{+6 zIwi`z{jufmmhm4%>&NX6qxEC&$5C#7aer^e)z|mx$6K!Ld0gMK)qykJH8*U-rw6^{q=a;U-#{{Jg&1pG?Rnd`w#F z4r`wk`BBy$6ZMbRzdc&tkHfcjTjZ~g+T*MpAK#UhUv1T2m#^=~7w^y4^Zj2N?a#fh z@2|&m?|+`H=eck1Qp>%cx4pk_w_NppzFxc?-;eLt`+Gmn{kVNUp8NXV-`lbG^XBXO zdG+(`=iS%y{f)5K@9&#$-`Dqc{J#7C-1~Lo`}6I4+pjx6KKH)<9d;bv-@PB7U-#}u zr~Pwd{lIN5>WS^R{o_kGoTJyuMxcaev>AuNT+2zkAQUdLHk`^SJtclQV=h4q=xqiM~_wn}qeRc2eyRYZ{z3tbH_xIe}zP=x?UzhH^zt?#CemuUu_v&9@x5Mb)?nNR^N}>lcV;esQpdU{vv8mjoMSJ?fah|wWmdG-@f%`cuf83xqinsA zRy~i~p8I)qAFuE4ODSw`j`uBkUt3SJ{aXLl)1$UuH=`|&udDdF^X+^8Nw%K1$6kF@IjMp?TQt>@mWZ^zqSJug?ipI7haZLgk}>^Me6$2rp4 zlcMAF_2TWs=OsRGzMj|ky!-ny&Gt7U)bxEEVtM3Sdpo`kM_T`qRllyhf4sl=y7cST z*Z1q%*Dt%jqTPpCEk<=cYj?J~Y2-JL+MS|y`=}jnFK+vOympP&kN5BWymqnv-tG{! z*NwL4{d_-hfA?Pf^~G&pKdx~<-@b3htM{AN>g#!px9i)9kH_1-p7)E}o_qD}`}&^y z_2TE%bMNo%_iy&D{XF;mc=bH)=ebw+arN!^dG-GOI?8$P@8`w0>*uYZ z^XlvS_WgZ|-#33>ef#cxeZMYzfAMwW*OhPIuRGtKUzfiA@zHfq?z;B%N2dM5_sNrr z_Pugc-?H{sR$njE{=)J%t$N$npJw@Ft3OzM*XlP`{WyI6*uQ7}{CIqS-j4kg>+jop zKU)6-Yky|-!!rGS{n$^9*8kYro_qDz@B8z9o_pJ?ujkcsukm`mUC(_zufCqYPS4A| zE-+w^UG6&e)#e__l`2{a{X_&{=VLw*7oY&+i~^%d+yuu{dw-y zbMNOh-p}239G-i>36}eI;_CbH_1ydVal2gZ{d_yVUC-m{$K~7e_55}C`ks3~uf9F+ z=hbuX=XreIJ@@MS@wQjrj(b0D-(I{w_r9LDAukYUT`1swwZQJ+P z;r+dO?)&xJ*Z16O?7g4o-p^|}@B8)jyn61(?dyB>*XjK|_w(Ypx83{k#OK-Hhq%9Q z$Jcl7@3VW~o_pV(zYafM-;TF^eQ(F!+rB^V=k3`0{-+l0?<}njusSfx54U#TGJcTd zzP{IfQUCr?zlOZ;*SF{EA7OSu8+jn2Czpv;0y`SfCyWIA?UqdzCPrN<%UVVGz z`nhkY#_M~#-0Sr1`EhtZufD$b^S1Z*cDeq(p1-f<`n!*>OYh&1_kMnz`tkYp{JQqz zjJ>x96z%W$t&WS@H(C2;tKRm0-gZAW>R+xMukUT|@742o`|jiF{l=7e9e%tGRo|{} zC+_dwkI%jD-`ie&J2%>KdE2-1U#o8?-mdT0YkXXO+`eD$=he4c&U?Rd)z|aw`2IaF zSG`|+eBR%4->wDX`<8@py{P$wht^N|VJ^#`2pRE2|rl0$GJKo>D zx4n9QuW>)$zxVg*dEDQR+jFnpKlTl^eSPoe_1EaU`0=~<8XvE>z54dO?dy5w5jq`GC-@mWt)%*K#dG7r^_ZsKk-;dv`_xI}kyzTq-{_Z{Z z>iymO{@us*m!kdm_E8;cZLi}ZAGdEd?=|+`F6VteUcLYL2KMXy;_B;pjlJ*3Yq@^6 z%;W3we(`?d*BiIv{@(WE^6GiH8gDml`+DV$)4i`(uAiSLKQ6DnzE|(()$@6+@p`@= zukm*LxZ~X0?!9_muKIqw*7dXh-$k*p>VNx@YdjChX=pqTWyZrC|Nn=r!}%(CjsO2K zo_pQO`q!_8(w4Q!Kcy|rd)=a5o-b`-KCV5iUk~$f?QXuid9RyC`R3;9*L*3o^QFzq zm#and621zoz8vuIB64Ldn}*%$KX3&DXDmlD9W7U#@mCU%ytAyuGe@uO03D zg?6FT!F*iXTfeo;dtK9Bzt`2w&%IWZRx)3|7A&tQwKMOvZL<9trBG*<+L*6jt1YiA zwK6~VT3HG@SFV~aTZ(4NR@-l}eyvQ~cHG7KwbI(v#dcu_X|>l~wDVB3^A^<#YZr=k z-lJNPw(Wf@dbRhh=(U6G->bdvMQ_`6Q>u zF+cd{zMYqnf9`t?=h;8sbJ6E}zSON&|J0<_Kd1e3+dr2#PWsiBHmX-^($8JH9!v4( zZn)n4^EF=zAbF%*DVy@KA&P&_mxMq}=Nb0P#woKa7OQFsx zh3@|OR&C9+6vI-r)m+imvpT(K+fS?QH~l<{KY!xS8~jgr;kS;hd%DS)*ojD8%Ozc zeAORM1s^9~8yGM7QW!5@bJ2Jy*m$YCuk!Jd*6@a>GPt!kLf(<(57dn|Hj<>v+jEFu&wV{v2Y*lg=|fPHddpILZIV@ls>s#m7x;I!=6^@%hE? z|2GWhJDq3Z@#5noRG(jJZGOq8^Gdj0YHgg}kd!sAx$vdML`e(8O$Z6}>qY(F+m%8wTxC-yqZjF;MU+=TfhAKj1o{E{m> zUhI8H?<;+N@$up{jFZ};*YHAX{^NMbr}uSpA1}G|zRl+un-9v&FLQkFQD%Pe@lrRh zgzeVdhlTz5{KEK2-sgNN9Vg)pt{W%Tp8Nb#H_zAwlFlz-UdgB9#>a`@kA-oPPw&g@ z^R0pLl27LszYnwPB7MFM<7IZ4@sf+?mt1K@dwp$_S~tJ=cnLL^&MP)gObn|{w>l$? zmb!m;3x0;x>FK;;>!YNSgp-7dAkpb=a+K#W4WR~7p3=Q;e@5*#^)9LTuezEr_eQxpR7<*r% z`6V7Nx%9rQ{{7fozF!C@G#W3tqFtx+JilyZ+lik48<<~m_WV*;ZCu)OOx<%!{qKJ> z^Gx`xtea=7?ek2yFZKDw$4lLFzmF3iH$Ja~c_v@-apE<+Pq{Mli;okpbB>pAz06GS z%Y2-8_2(9^;d;!KytcFP+BP}gb2_n zpZDvk9Y;F9%>DT#SMuLQe4Kcl`*_Kv_hUXzeBAi`*uNVuxsu&Cr1OeD$JqD8>GMm6 zaK6&}vUq;+apE;xhd#d)O6`*MD@tvX_sM@BsBB<%ND#<19N~e15UlKgaxH=gIEFa%Jwva`X7SpHJr*A16L;?1dKXg{SjM zx$zR;kNJ2h_x$4bVP1Wlc=d4->c4&7_vaU{b>k%LuWo+vabw3@_Z$=6@G|pD-S=WX zZtOgjy!vr{l!Ojn6CLda6yI zTk4M&pI`ib%%5L;yi}#n`*xlCJToyI-VCqdY}Ng{n|ZG`zZBAO>kye$3wY z2F8p3-qX|Zf@5BCY$4eod zUwmHiT7Q1YrO)|3ZZ-@TSad%Y&o6cNWAS*Ye?R8mkLB#SzwSA|+<0l1ygxo}e7wZ- z%XGgVD}}c+H_^t0HPe!g6HC(ju`uqYH!#2WIPp5?c$xeCn14T3cRyz5DH<<+A6EB0 zncWW;?R`jpcbWTm$)(>X`#fVW*v3nt%>452$BW;W`Ml!Q=ao?B{5^Bt?#J@;I9~kw zvATKY-_0-n+~V_0-sY7;nfayeIp4=kcpqllb4=a5;`5B(hlN@MOlxWLormT!mN6w~j;tnK$>UhB^<{@mifzr^!PJYMqtJEhk-p7*!1>%iw3 z8-LOKQun>G&nsU4-^PnS@7L9E9{-u+`!Rp+x97<8x!*qD()+RSJ{IhL%;%L{ss8gz z-S7K;9~R~rpI`jc!(%YG=G`~zq=a*Jtzv=hMZIZT+7oT5d+5D1AYnWf+@ly9Y zjLok$US`>NiN8z2=3ids6O}ufyzlc1<0BlGJ?Gbr7au48?S8EMcqxYKAsr`S-md#z zEZqP4@A@^-c!__1@#p<~>3{$ISp4@eHTHX@|8D8uk5$;;!_=m~>xXgU-zWRuk;V68 zKEL?$eqmxbt$sf?ExhWwe|HOhn$_txUi|wp|9i`LewkHl7q(&Jq)_zV`Bx1awEM9N z`@8MhwAPGO+^llu2#K5lB$ z@5TH%rvCg=_uNwd`Ng02bLqTd*G+UkHs}0e*I}9Y#s5CJCi*_v|33NmiQ)47#%e?J z%Qto(_PbRdFaF$LW6v*sKUO!t_;_h(e(~|*^NjtyQ1tg>KEL?$OI`JG;?=(w3-b&nuzU{r#BFGd{0)jmL}MhlT3n#s4m5Ugwv2 z{k^j8d$DjljEm#~sA=%`m(9Y3l0LWC`6$`vN?Lth2{o7g?h?j}&ntdE=J#Q_l7FvU z|9RiXOZ@lzb>k&mFSY6K`{6!pR{D1_;q$pB`g=?N{l)*D-=AaZ{ysVUPBPo?!>krc zx$uwGLf2`b{kYX)XcnyI()p$S{n+gEzN~F>UjHt(OIGyfmb&{fJO1?VVeCHB=9gJ^ zALhUBmwSHk@5kc5=P#t+Bl~#qabow)>3vxJ`Nh9ij>n7t9ZY82gzF}Kj;TLh{C9ok zouprl{qDklACsJqy5|-jH}%I$-8c!)@ALY5WkchoBK`fvj-%w&#&rYZrJ`h?6J^HB zoZlzM}$JlR2MX{5d9EZ?n?pm~bEF&-;avU9VREvG0@p zdB4`aANw<{;cfK$vRRK$3PD$^-NNNDt7!dfUimYvVI0;*-;e!MwBtzc$Lu)mJ}hVB zWtM#}mP@P8E8%m)pI`j@Vrj3`n=Y6{`N#_;6 zAM^ReYnWf+=NEsD$=NumtF7#Hd$plavKp771!=YQY~18*oYd8JVTEXZnfLKxuX~Q) zE9=iM^ZvZ=@mhbp#PiD>&-=wi z!|{}&>bci)-t)L$W4r!&Z9E^1=R;?palC} zW4rP5RcBy+zdqt|Fu$LV#@8P{7vk%pu|4%?oEH?|u;Uv&oN_v<4b2lM;+Xng(Qb0NMi z8ryR}FOApB+|Ntn^%Bm@{C<7J*TMXLJ{n(tw%fkxL4&+t@{ol}($*l~lfHJvNLI2k@Wpb;-8lX&WP0;^kO_|)QR0DU(|dI)O5|^+EsFH)V1gg8r{?%H%cz{a?wH$!!e!zmh4F+XVD~B~vEX z4T_*lZZp^%iXhh=dVn%H{U68le;HFIw-sy+%H+0zZ9$pbcCbAtliLAy1Z8qN!Ooye zZWq`Ul*#P|yMr>hJz!5zCbt*t4a(&9fqg-l+D3dz?4g_U#J>eiwCf5rNhJ!$^ zHyi@Wov>lRFNM2W4_6z=@zt?j#rq%H&RlQ(z>>oeHOcGP%>?3{WO_CY%Mz z}a%ZX%RG zncP$GG?YN@8F&_y$vp?pgEF}n;6+d-_Y%Ae%H&>wSD^`1K~t!P*Wh({115noxdosF zD3f~=-U4NEEuj@ClY1K`gEG0+uqY^#dk5YHWpYbE8&D?q9!vpca!bQ9piJ(4_yCm2 zwS{(|OzuON3hh8{MOX=x$$bPLgEF~QVKq=D_X$h`WpZo6S}+adK84RfnOuA50LtV( zhc7^x+6wv5+sYJOB=d!7wCbdASOx3|V%*P?k8xU0&mEI1GcMVMNICa!sK+WZC&b zS>l>eo8ym%0ZJ}Mrvh#(q#4S%<0e=VF4tK+Zkmco8 zgq1>;oiCInZe{8!`1|2LcnBU2Szc~cSS@7P`9fLZR;R9ke;giziBJkzUT#fTD`eUE zLRsS0rmlm37M_6@;Kh*T<=R7skY(o!Wr^!ZT^Iivybf=`B(n9O6RZy#ge*H>C`;Ui z)Xw;~;7xcNCWkC9w-Iy+S$4iqmbi_nUGeY2JMbP%30Yol6WBCl+4(|Q;<{0b_z&QH z_zh}1r;7j-#z6n`gZU@*gWZC&bS>kr0?u`EqzJ>4Mhmhsvc7a_( zmYpw@C2lwB?)abJNB9}0hb%9*2kaTL?0lgtaeGns#{UYxz;Ez-$ntXgz`h~N&KJrO zw;y$X{Gad#`~@>YmX|vK4h&g#zEGC9p45ZzGvROe2WEvVFV_nW4q0}-P?osf)I;z! zFdJ%VQOs+5Vg96aW=Zhqh>ko%PrMV{5D$C>!hXIz!HHB)+et|hdxOzucH$}+iyVG+yZj)oDI z$+d<>Et5M2j^nl*ye7 z=YcZ0Rbe$yCU-tu0LtW6hc!T%+=Va-l*z3LYk@Mki{N5VCbu@M1IpwsflEP|TzlvM z%H%GC%R!l3M_3n>$z1_g!g|mN)`tz?D!3Z1fonmT+=kE@l*wHO*Ml;-ji3uCle+;% zgEF~|p(`kpyAf^zWpbOqrl3r242%V3a^0W^%H(c_aiC0YGuRxI$=w3uL77~4=mE;) zZiU-GncNnzB`A}-9qs^Sa$CXHpiJ&gxC@lYZ3EkaGP%290w|N)4z>qna`(W!piFKD z*b$V;-3RxBGP#{#XHX{h06YlFhc7^x+(|GJl*xSwUx6~Yli?IlCigXb1Ipx1h0{Qp+_&%@D3d!K z&H!a{-@^}}OzuoL3zW(I2tR=`xwGLMP$u^?Ob2Ch=fZiQOzs!>6_m-H4;O$kx!>S- zP$qXFi~?nHf54xhOztAM7?jEV1v5aI+$C@+D3kjeW`Z)g%iwZQCif4_0%dYnz?Gm( zZZ_1wRd6+21J^=r_&KDIrxq-eyAG}oS>h^Kt_)drzED<>YXVjH(QpIY2sec+aZOpS z4q0}-P*#v@2F>wfVGP_1<3g6W1z2tovg~}JtRS}_EQB8qx4^A%TgVdElI2z*%gz_d z3UUj>BKSMtcDNJn3R&V>v%F}?vh#(qg4|-TIDP`$4fnvkAxqp6EVl_+cD_(nkXsU# z!ru?~!2|GM$P%|S%gcl;J6|X($Sn)Y;U9*F;1PHRYS$4iq zR*>5SHpNeY_uzf_AY_T_#&R)a+4(|QL2fhH96uF4gpc6kkR`4=%RNGtoiCIXm>0dFzm6=L=;8xou!u{1@;!dB^M$g4+zzlK{u}rj zzJ>2Xmbjf*-Z^C1`9fJiZWq`U{{wsvKf+HTOWbZO?;f)3e4(r$w+HNrpAJ96FYs%~ z61Nx2dxtDLUnncc?F0Mbe}~`T5BM`=iQA9m{X>?WFO(JJ4uAvkGvF`y8)k+qaXnc+ zC}i3BLRmqs7aWYA1^>Wos0mr(db4~8SuHgOd0?T?2M&e4&<~W!RX`;ulj{$Mfik%! zPzB254u=7tOs*+ZgEF~+FbI^%HG}4$Ol~j?0cCOvKnqYNcLWRtWpWF`LZD1;7z_tx zaxI}1D3d!9jsj(J3&SFyOzvnH0m|fB!=j)}?ie^0l*ugyi-R(`Dpl*zS)cA!k|bT|W)$t@2n zfHJu=;Ve)lw<4?r%H+<5b3mEg%CHJ3lRFpA17&in!fK#Q?tHial*z3QYk)Gj3t%!Jtg;1$YsZ z$@PXqK$+Z2@G|s)L!mG9gIC~Hcnw|$Wpe%DFi24!-CVF)Obdk5YHWpYQrP*5iK9!vpca>HOaD3g01J^*EMN5WB{OzuON3d-b; zh7q7l?j!gZl*t_f$AU7sPhc7-lRFNM2W4`f!e^jN?gTgyl*xS#Uw|^XlVBt$llv0B z0%dY1!zrLl?rZo4l*ye6r-3rLZ{a&oCU-iV0m|gQhaW(h+?j9|D3kjUegb83XTv$5 zOzvlx4$9=th4Vm}+%NDeD3d!ME&ye6zrpXIOzuJ$1)|@M0Y-^Kg;vvh<7s?X11hozRUbqMDgZo2PA-5zf6|(Gnp)7GrQjw>+#6vg~}JEO9GRSHeFDPryVdg{(qu zWmqL-+4(|Q;#Q@uhJPBKf@k2_kX6X74r_!gJ6|YE+?v$2@Xy0@@B+LTvI@DiVV#g= z=L==p3@_e4#9H z8&SL9C&SzD4!j$(3b~D;Ysj+mg|ftLLfsTU1>S@A;e(J>$aRBa$g=Z=vczph-5fs^ zK7^0pJ@{?0lgtaobS0#eV^x!oiCInZU^d)_;28A_!hnkS%utAuye?=^M$g+?LyrZ{{wsvKf+HTtB~6b zb`M#0zEGC9J*a!)r^CIGPxt?%P$g=Z=vc&bG9*my_|G;dh30Z|)Z#X1m+4(|Q;`&ez zg&ciqsd*?^Cf67GStdtE%@;`~*B=hEOs)x3StfTl46sbDDO6h~HxLF{Cf5v_TP8Oc zhFB)I0JN}7?g$uancRY~kY#eiV7O&+EuocVa!0~ZmdPy)i&!RiG>ouJt~D%bncOjO ztYvbG!Qz(59S6r-CbtB%u}tm+IMFh>C1EMcgg?5(7oepPMCbvASV42*RaF%6qE5b^a$(;@7SSGhJtYVqmxp1Cka;w5>mdTwD z7g#2@I;>%t+=Vd8GPyNjEz9ICf{QJaTN~D~Ozsl6)H1pD(7`gf%iwb82Q%H(c@n?RY| zCa@_elN$qLL77}PD1tJ%n_(O%liLh72W4`%z<5w5*ByF*GPzsfHc%$F1#Ah*XU{_Ek_Ygb`%H(!~-9ee$Bk(9FliLIK1Z8rM!Q-GzZZFsyl*v5-Pl7VJ zePCZuCN~jEpiFK**dLV1Jq1sLGPwibKu{+43_J_UeF#%QncUGZ0+h*p1RsMkxntm1P$u^YOao804IVnxzFJX zP$qX0j09zJU&2?QOzvbj1(eBs4c~w=xl`dZP$u^+d9w?Lh1%3r(a_7SZpiJ&J_#KqVT?nH-S>jr<+$v<*`9fKx+`_O3{tmbu?u5HSmblg| zFB-D!e4(sTZZTLKKLPHBd*I%XC2k3p+k`ATUnr}TTN0MS-w*e}1MpzT61OzV%Y-aD zUnr}TTNakXKMW7SBk*X*64#dHb|K5o7s@K-mWLJakHcf|1UwnC#I4BkN+HY67s@K- zR)$sZC71|L!P6m2+^Q_E7P9Pop{!DFbyx%cEIb3x!Sf+Y+?p(}6|(Gnp{!DFZCD5Y zBD?@E!OI~_Tzi%~ge*H>D65p~2VARYS$4iqRw=g$Y>J-(@4@@k%Kd{D@;^M$fXxn6KEem2xV zE$tlihC`qa913KK%d=dtOs+5VvrMi6=$vG7{oydnLbYXb17VP5 za?PN*WpaaIh-GpMKnu&{j)0+-$t?&AStd6OhFd1r5?Wa%cO)ETncTv#h-Gp|!wAdd zTEn82$sGg7S|+y`EN+?Had5n4a!WuP%j8ag6D^Zl5|*+|?j#s#ncUK_jAe2s!zq@@ zEep$8CU+{F2Fm2xLOW0^nl*ye7 z=YcZ0Rbe$yCU-tu0LtW6hc!T%+=Va-l*z3LYk@Mki{N5VCbu@M1IpwsflEP|TzlvM z%H%GC%R!l3M_3n>$z1_gf-<@Fpc5#Qy9%xbWpeAo2B1vt8n_me$!!RoL7Cija6KrK z+X%XVGPxUIG$@nX7`lQoxf|goP$sttYzoTc#=uxmCf5y$piJ&&7zfJaHiOMUncOWf z9+b&-haR9z?pC-Bl*w%YTY@sV+u;sSCbt!A4a(&1gu6hQ+%~W+D3iMzCV(=z?O=OQ zCU+0q3(DknfE_`Z+b_Qi~55R+ z_JBP>ncQRWI4G0b3-$(Oa!hc7^x+(|GJl*xSwUx6~Yli?IlCigXb1Ipx1h0{Qp+_&%@D3d!K z&H!a{-@^}}OzuoL3zW(I2tR=`xwGLMP$u^?Ob2Ch=fZiQOzs!>6_m-H4;O$kx!>S- zP$qXFi~?nHf54xhOztAM7?jEV1v5aI+$C@+D3kjeW`Z)g%iwZQCif4_0%dYnz?Gm( zZZ_0_GP$eZYEULuo69xHStfT4Tx*$J9tt7L&KJrOS3#}BUk}&84KO-nP2|{kld6zq z=L=jq#7scNNcf#E;A!JSD7K6n@mYpw@C2k378~nX+ z58Ma$hpdU*lCV_Bvh#(q#4Sx-2LB*D01v^#A!{PHEG!qY?0lgtac!yX@Q=bH@EANE zvLlC_5p&=J-RS$4iqmbmq(o$#;2 zYw!k43Rx4m^V`^N?ld3uTGxPVIsJ7(RkeU|PtU$ZY{z zhAcZ@C`;T{)UENK!Kd&!d=auHa@)YRA;A!{Oc02~;y?0lgtaXqOA;b+3% z@DI!iSrfTlaB#@7^M$g+^`;(zuYuW6OFI{`CUSk?(2ylA&vL;sxxUcPGPw$LbYXb17VP5a?PN*WpaaIh-GpMKnu&{j)0+-$t?&AStd6O zhFd1r5?Wa%cO)ETncTv#h-Gp|!wAddTEn82$sGg7S|+y`EN+?Had5n4a!WuP%j8ag z6D^Zl5|*+|?j#s#ncUK_jAe2s!zq@@Eep$8CU+{FW|>@DXlI$+>2QW+a?8UCmdTw7 zXIUn5!y1;!T?nHrlUozkvP|wG zxY#newP79057;T?f~LGP#YQ3n-Jj0Y-x|xs9PKD3iMpZUSX;o4}@^Ol}N}1!Z#Gpa{z3ZiaE7 zOl~vS9F)o30^>oMTzBXJ%H(c^+d!Gz7O*8Kle-=60A+Gp!PcNm?oPN1l*w%a+k!H= zyI}$-liLoq2W4{iz`dYMZU@*Al*!!(_k%LIonU8BCieh52+HJkfn7nF+(Ym%D3jX_ zb_Zp0kHDj#Ol}X@6O_q429JX>xxHX-P$u^TJPFF=_JMstncPGufik)MV1H01_Y^!0 z%H$4!13{VGGw>`Zlj{iwfik)0;CWCc*9#5?WpXdTi=a%dHyi@WX5qu2F zeFr$OgIaa$^8gFfik(X;T%vV_cKffWpd}j zd7w=07x)#F$(;`ufHJw?;CE0acOi@dWpaPOpP)?cBDffo$^8X0K$+Yna49I0`x|D$ zWpFuM0awC5Fbifw4Jebl3a$oaa<#c!Rn9WGYv5YTW0;8wUD?g&|x+`_O($g=Z=vc$EfE{eYk?u5HxLddG*7K6n@mYpw@ zC2k378~nX+58Ma$hpb9&Nmwdm+4(|Q;+Cc^gMSbnfQR7WkX6Yo3(JKpJ6|YETw7{8 z{G;#)JO+=4tV(WqSRrKD`9fLZR-~?ke-fU6iBJkzmE6j(O31SFg|fu0N?i^AG&}{* zz_TH%l3N|t2w8T%P?oqgscYe%hv(o0crj#Ea%;mnAj>+HEIVH)OWb=x5j@4pTg(xMaZh;wt;O!mYpw@C2l+F_V};h zOZXbT30alg4zOd$vh#(q#O*}g8UG!83*W;JA*+(x1$GTtcD_)SxZS9`<9~u5;b)j0 zvMRYfV9$_c=L=>IM|e4#9H`%(AD{|SGC7^r@xhpnHmdULR>sTgt30!KKTzlwXncQV?xn*)4VO`7Qu7E2olUoluStfTCTy2@$`mljz za@WAMmdR}hoh_5Q4z9OMZX@VoncNL9+A_I~p{r$bH^NPp$!!9gS|&FJ##$!V4T_e@ z-3;R_liLh7w@mI97;l+ecj#f6+^uk%WpZ1z_ymj z-3=2gliLoqw@mIHxYshd9biYxSKw7pCf67Gfik(*;B`6SXx%1%yP$u^q z{0_?GE`(8_Ozsc(6O_qa1Q&xcxxZirD3iMcE(K+Bf5S{rCU+TJ4$9>IfmxtT?h3dP zl*!G88c-&86z zR%CglkR@};M3ze-Ybv)gtP-+hE_sUOr$d&wRasswWXW9e49m}ktf}1Uutvy|x#T&P zpAT8$)?|4tvKQb*cnMx6TN~Da_Rt|@+4(|Q;yO~-#lH%#z-#b&$ePNn2c1HeoiCIn zZhh(o_(|{vya{iGtf}0F&^ct;`9fLZHllXHPlmVQ9e6inP31O*t|80L7s?X133XHa z6nGEbhYv#5RIVEoLzbN{lqGI6>gM>V@F9EzABU`|TzBXZvg~}JEOA>_FO(&& z7xiHLEcgdzLrut<%JqgrLYAE`lqIeY^-##sravYAj{;ML37LG2E!1`GGP%>?49nz}hZQW7I}^^bOm0P3 z$uhaK;T+54R)$q9lRFpAvrKMPSj{rI^Wg%^Nf3OzsM}5|qiU2c1Bf+*NQjD3egHwMOnGP!P01Z8qJ!#GeTw;5~> z%H(c=@t{nuJM;i$a<{^5piFKH*bmJSdav1qXvNxfkF?P$t(K4gqCyFTu;8Os)?c3d-bOfmcD9 zTwmx1%H&>y*Fl+Fe>e=3$-M!SK$+a(FaVUvy$NrDGP!{;2$ada4U<8c++Y|2%H-aG zcR`ul5ik^#$-M_tK$+Yy7!Jzh-iHrBncR_Z6eyGX5T=4MxuanOD3kjLJ_coS$H1|m zOzsny2Fm1)gX2M&+^6svD3dz@P6TCgpTifROztEY3CiTYgs(uE+{thXD3kjdz5!)& zr^0EVOzvCw4wT8A4rhQex$ofzP$qXKoCV6{euSSuncUfM4k(lR8K#3WxpUz>P$u^a z{0hqC&W8&?ncQ#iJ1CR85JrJAxj*1fP$qW~Tnx(O{(>2B30w-7!R7Eb%!Gen7ATXu z0c_@S|J6|YETm`ife?43WH^AtSRm(Mj zs*q*p3uTFGO636J`9fK>Tr+5n9}8pPW*8T;#4W&bi;!jK3uV=E3&KM9@o)>=3b%zU zaV=SH6|(Gnp{!bNVORuz2iy*K!d)RtTx*sW4Ow=+P*yFs7%YyT0C&SZaBs*Gw*<>= zLYAE`lvT?u2}|Mchx_0Gcrav%Tbku%LYAE`lvT?u3(MgjhKJx0cr;{*Ys+%GkY(o! zWz}-a!wUGv;W2muo(x&yR%CglkY(o!Wz}*k!z%a^OoXT4>5wIERhCx^S$4iqRxP(W ztbu0MUVxY2<&Y(=J5J6|ZPmfHk2#ZQ6v;C=WYWQps>axrAt`9fK>+-9&jekyzjAHl~VOI&xB zdxR`IUnr}V+XA-4PlHe3Q}`@oiQ9_htwWZbFO*fwZ3Ellzktu-OZX~eiQA6l?L(HG zFO*fw?EpLCzk#peTlg+yiQ9?gokNzLFO*fw?E<^ve}M1dNBAjZiQA3k-9wh0FO*fw z?E!n@r^C9w;(KJ zncOfKZkb$5Xl0q)k#Lk{atp&EmdPCrBP^3^4U1YPcMKeBncQNqxMgz3!SR;KEdgyT zlRE)Uv`lVESjsZFlVGG}a!bQ9mdTwAr&uPpEG%c4+^KMyWpZtyon>;T!x@&zEe|VL zCU+*BWtrTHu##nRXTv#`$*l~lSSEKaoM)Nbs<4`6a_7SZmdULSYgi_CA&jz2ZcSLr zGP#T3V$0;#hIK5Hy96$^Os+k2uuSeUxZEjp*3 zmdR}fn_DJ#3yil+t~>OwOzu{=%`&+yU`xy7ZihQ8liLcmwoL9$xXUuRZD3o= zEEPEwMQNOJDsUhd<59^gS9;$a@) zQ6A%QQt$*RNyU>q#nU{)vpmQ1yugdR#LM_MEP56HRz_(sCF&ce|$71Ah6{Y7JGT?C%c??A!50S@8WZg#| z6Or{DS^JT799iR$^%_~rk#!kavr$g0&B!{7!Z@v`$Xbi6o5-4qtdGdriL8Uj8i}lh z$XbYsQJfN#q!gtoLs`mEo(fc?5|yb!RjN^)8q}l~wW&j0>QSEtG^7!YX+l$)(VP~v zq!q1cLtEO>o(^=R6P@WoSGv)i9`vLaz3D?=`q7^O3}g_48NyJ8F`N;MWE7(r!&t^K zo(cTIM1JKre&-MVf zMJ{ofD_rFo*SWz>ZgHDC+>N4`#QZ0gKoa64B^ik%=N|6mKJMoM9^@e&<`Ev{F&-xc zPmq#SJjqi$%`-g9b3D%ryvR$u%qzUgYrM`IyvbX<%{#oyd!!}}X-UWXe87i%#K(NX zr+miee8HD|#n+_g8#3@M-|;;^kdaJeCJR}~Ms{-WBR`RopUFjT@{pH&YE-8NHK|2y>QI+@)TaRrX+&e1(3EC0rv)u( zMQhs7mUgtK10Cr^XS&dpZgi&yJ?TYn`p}nt^k)DA8N^_QFqB~oX9Ob|#c0MbmT`<{ z0>3bkU-^yS`GY_Ci%Cpo3R9WJbY?J-EM^HyS;lf!u##1*W({ju z$9gufkxgu73tQR7c6P9nUF>ELd)dd|?B^d2aF9bB<_JeQ#&J$?l2e@K3}-pVc`k5~ zOI+p(SGmS@Zg7)Z+~y8<|MRq;5Yq&b5GN_gNF+J;a4+|9KM(LA5AiUM@F84j-r{ZE;a%P%HEBpoI^O34KI9`l<`X{U zGd|}FzT_*uCOzMffp7Va@A-j@WFj+J$VxV{lY<}miJbgQE^?EHyyPQ41t>@%3R8rl z6r(sLC`l}c1 zg^B#iZ~V?5{K;QTVlq>h$~2}kgPF`?HglNEJm#~2g)Cw*OIXS>ma~GDtYS55Sj#%r zvw@9lVl!LV$~LyMgPrVRH+$I2KK^Du|8RhV9O5uXILa}ObApqc;xuPC%Q?<-fs0(? zGFQ0DHLi1mo800ycewkXm;6a$nm`ibBqbS%B3%tlnyv!@S%4@vN8@$O|yv;kj%X_3I4QWZo`+UHMe8k6m!l!)3 z=X}AJe8tzK=NmHcE#L7yKai13WF`w)$wqc^@FPEwlb^{&Zt{?qeB`G91t~;ficpkd z6sH6wDMe|@P?mC(rveqJL}jW_m1+= z(3W!$9XPrkxN|W z3Rk(tb#8EzTioUjcmMNtKOWNrk`N~;$w(wQ_i!)waX%06AP?~{kMJmu@i-}Xf|R7< zNuG+l;K0v`{G700k*T zVTw?cViczYB`HN|%21Yal&1m}sYGR}P?c&_rv^2tMQ!R(mwMEv0S#$HW17&EW;CY- zEont-+R&DEw5J0d=|pF`(3Ng)oEPH>V_oaPK?ImdY} zaFI(~<_cH2#&vFRlUv;84tM|a(Es1(49-IazU4c<=La&9iOggnE7{0S4u0e(a`H2| z$W0#dl8^ippdf`POc9DwjN+7_B&8@#8OlHNAm8eV=s#1;W)SxD{s7)Q}QjhvH zpdpQDOcR>YjOMhUC9P;p8`{#2_H>{lo#;##y3&pA^q?ob=uIE`(vSWOU?77S%n*h$ zjNy!6B%>J37{)S=@l4|rna_?!Lw!vPL*h{GJ=D91R? z2~Ki~)12Wf=Qz&=E^>*>T;VF$xXul3a*NyC;qHH)_W$=egY%GqZ~2bz`GJgNA~RXY zN;a~SgCF^cocv5Ka+8O=lxi$tXrMhOvxeJQMhZiTuiM{LUZz$zM!jGEEMhTBSjsY%vx1eZVl``6%R1JxfsJfpGh5ioHny{ao$O*ad)Ui9{$@Y_ zaDamx;xI=z$}x^}f|H!$G-o)=InHx|i(KL|SGdYGu5*K%+~PKOxci@%{Qvu$!FkBQ zw|vL<{6I!Bk(n%HB^%kv!H@hzPJSjAxyeIb@{ykc6r>P^DMC?-QJfN#q!gtoLs`mE zo(fc?5|yb!RjN^)8q}l~wW&j0>QSEtG^7!YX+l$)(VP~vq!q1cLtEO>o(^=R6P@Wo zSGv)i9`vLaz3D?=`q7^O3}g_48NyJ8F`N;MWE7(r!&t^Ko(cTIM1JKre&-MVfMJ{ofD_rFo*SWz>ZgHDC z-2Knn{di0hNJ5;XBqNdJ+{3-x$NfCOgFM8;Ji?m2~v`ZCwVIRzxUDczAl^L z|9?LkH`nCOc|0BF8s4jhJs9!#G2`(5HSE`je}H+1_q1UTNBl#~KD^IOzQE;=unxj| z->?TH{xOpr*LG~y@lTlSxRzt{kAKP($2A>W6YbQnu>m>d;(;QcKJj>-TnC`fm zhdnUyukrW7`}MG2CY~POH@t@rduZa{;OB<-^@%#X86N9glQ5 zH|rheb3DT3JZx~B*YR+dHNMfY>%$(fxIY))FqR*m(_m=;F?Ee2-)2J?yuOn?u}rk3P`7YF^>9nz(Zv_T+`f zXxO(G9;aasUwEvBeSYEb8ut8!$86Xa7#_D_k6?K0hJA$L@f-FOhR1N&cNiYWVGm+> zEQfuH;qe^yEQZH)*w+{y*I|!icx;D#km2zi_C$uqc-S`?9_L{XWq7P7Z?2DPGfkmQ z^$BgFDK*^wRN6?N)Q0+$HqfWFzCNS%^jWQ|&uJZfUTf!qviBnEvxTo8BMLFHI0_iwCppJjD-6)-t}e{ z{{0X8hQo6}^5VKl%|1LAgneb<`5^2u3(pB*A6j@`2z%1Pb3@p-7M>r%9=7lt5%#%- z=ZUcAEj(9*eR1LWBJ7b1&lzDKU3lIId+Ne-N7#3lWWVOrfAl9kpg-zC&7p@hyB^kT zdPKA8QO%;qG_xMpOnO4y$M9Se_9cero3KYQJm-XcjNy4F>}d?oJz?Ktc>W1{Aj5M| z^7MK^zt)TTm0r>>^|F4USM+nes-NjK{Zy~(CwfCa)|>i~-qH{Cwtk>@^nJan>CC}6 z;?70nwua9%_hSwQk&ue(C5%$G~*BfDvY*T3{Routp|WPL`b=+in?pVDdiq)yjVIzvk_?Bm+HN`Oz+X!Gm6J-j9g z`{2XtqVWHd3a^dg`=coOo@OqG^Xg)F-I{n;tEpLr*RP3Y8(zaEnss;`n`rjowQQob z5MIwFS{vatZKAajUe_jCJK?o$;uWp0m$iXj(uR6b8|ejYtmn0fo>OZi@vOQ{iD$IA zp4Jw6N?YnlZKWr)wI0_tdQ98uQEjJ3w7nkI4thvC>Ot+K2eh;Pqg`~rcGbVNo9@%@ zx>tMX9_^{SwU_SF-nvu!=nn0x+qIu=)Bd_u2j~_ZsGD_=ZqmWJQHSUT9cmVCPxKwl zT^{3lvr6)%nq88w)GU*Ht!A4fy_$9SU2LM+hxbW|)DAPF48QzP_ya+&8G7;yUx=bI#++xIr@{%)|@&^f7Y3rOJ`_qovwLwn&#E1 znop-_ex0lZbdna-zqFA4sfG0qEuz0`QT91N`Cu#}(MN8@gEv4hNw2spYnS%u&0Cz@UO-F~83hTrce znr-+Uf1+83-}5J$efV8}qO}ly-=Ao0gx~omS}WoA{)yI3k`&ropU_^KQhRDD?V(R< zcYR8`>C@U(pV2P*tajGtw39xs9rXq6pf75BeM#Ht%i307(Kh<3w$|6QmA8*3}QRj((`M^&_pNA8Sqh zL~H1$T3tWWYWlfW)i2EApqXUA?TwFey;&vfRB$YPisLvqXqP= z=GSwYPtR*!y`XvYqUP31noBS1&w533>Q((ouj!9^U32IS&8{~!o8Ho_dRw#T9nGwF zHIwt_v?XTLnEs#%`n~!>;&&R?Z#AiA&}90JCTe<3u3zgt`jy_RU+R5k;oL?3Y}Cr- zk*+tZ@Ht_k*@e#v6U{PwUYKaM;WNWTvkspdCYpWt>@d+<2%jG&S{var#6)W)e2$oC z?S#(~6VvGwO{-HijZV|lI$ht>8Tziy)OU22zOA$MEuEup>Rf$8=jrP@UtiM&`l>F} zS9FoStc&#}U7|1QQhh;}>GQf=pVJlktgh5&bd^4>tMw^eqfhEuO{MEJrLNZ}bc3eQ zjrzE5(#LeOKB`;v5#6c}>o$Ezx9fwtLm$wcdcW>6i^FD;p4KjpaJ^Y2G*GikXsBkH z&`8ZTp|P5ELK8Lngr=HQn`vB|YZ7gt3EEO)YUci0l;CG1yQ^*Vj<(g?+D>n2d%dY0 z^oDlS>)J`LX=lBvUG$1})yvvVFKKtZs6F(8_SEy*OV4R>0mvoL-dFaRp&2Zm>$yMdQhFmWCwJl{-dLGzmC?wb&T%QvAS2s=^h=g zyLEzD`2Nw?I9KuTea$K%lbT&ZW;M%%ENZq1S=Fo)vZ>i8WLIk;A%|KU;ddX&td;Ql zk7U+PLQdVFKkItUrRy}euGKucM)T@w&8Mq0zpm5*xS8UX zi?p~d)DpTtOX_?rrSr73&ebwHN6YGLEvK`zyw21LIzub!bgiV*w6adsDmq20>SV2^ zleD`2r8V?Vt*L)#E&W|<>u*{|f7QA=QS0e1THh>=Mp5(?ZCxJbdb3JMsb-guO3gCi zNgbKiAg!g|^Z!wWT_b$y(^w+Fa9XGyO)JY6fkh-)dw1P8;d> z+E9Pc2Aa_$cCn*smN9=Onc2pUt69fRsM*I(sHS5?+HT&2swH9KtHNDPJYb7>UztVa7rOsDBBep<4*M<6-F49kR zv3{aU^kZGBAL%mvP?zfmx^UI*0uVcuG4pQy}qLx^lja!Z|Np| zQ#b1yx!j97tg~7>u`YVQcGdf|o8GJ4^&ah^$+f2@YA;Quy)~)!(YW^2B-&3Cw7W>V71>twZ#d4%M4FOmFCLy{;qlnvT?~I!dqTXuYgs^pcL%i#kp( z=y*M^6ZD+^qGxrYp3z_RwEm{2^mjd}f2i{v`%{nWUwTX@=~11mM|6r=`18?cIJfa| zzsxF@OAl#oHOp8YHQQKTHS1VD-LLu8T8I@;Ya>=rt(90I-J^wdw-(V|T2yyxG2NlX zb-R|(ZCX;dYAM~KrFFBG(M?)bH)=WEpyhSFR?u}?QP*lEU89wCwN}wpT2)tSHC>_A zb-C8iWm;30YAs!&wRN%9(M4KU7iv9Sp!IdWHqd$6Q0HnRouiF)wl>jO+Eiz1Go7K$ zb-K3DY1&e!YAdra&*)P+xjfMIW)*u$%`Wz`nq}-2HQU&$YSyvW)a+xgtF;h&L#>V2 zoBFH1r4#jS{YBr=3Hq*%*Y|Xsrq;2VM#pGc9j)nfl)kSc^#dKDAL?-ZNQddiI#fT= zA?ki49ju?}ApKn3r`Q)dK)=-f`jz(6ueGnH*FO4<_SOvAOTX2g`knUB@3p)Bpxrd1 zcGXPUMKfz>&7#g{EUR|ZY}!GyYkSS1?es@&t3PQQ&8e;RXR|mRMUndz`F?S?J!vbK z&#T!*7t}1Hi`rZ-saZ#t)$F4yYAr-p)!K-zskIVa*M@pS8|Y1~ueY?G-qyN$N9*Wa zt?l->n6#FIXbnxG)itiwG^tkAWLiZNwX!DHN_vl0)O)pp-lygDel4dDXjy$w z%jiQ|S|8R@`iPd)N411Lrp5JfEv6~7s6L@ZG^G~SR9Z-%)Pnky7SN|Pzdob+^jXcT z&uJcgUUTaUnoD0ai!)IaeM}dZ`@7z(q6KPp(Ly!LXpx$2v{=nLTB6x?sagxsGPO3M z2v2>T*BVn^n|J zU(oJqmQfEi+o-3Sb<|7EKI*O3LexjCji|3$D^d79sWi+vr4G<1bfBitLHf83*2i>+ zKB`0Y5gn!v>u`NYN9cn(QXkM!dcTg=`*e)nt7G*Z9jD24ye8@dO{TwSQk|%A{Z*6b zZ?&@E9M1Ca*mvoL^ z)VX>==jnN!uVxc{h}#n{O&^!t*0^&Xx$SW?2sbirCXxFVH>1dXjhk8IzQ@flaxUWL z8aX%dF6vyxJF9aS@1)LUyrVj|@eb--$J?uOA8)5-5pSzz6K|tt6>qI(7jLC0w56JD zyoH)|yt$fvyqQ`H@uq5R#G9zK5^t=}Xd``A8|rh~K%dw8`hwQe7qzZhtMNMevewpD zw3fcAHT5;Ep|5LoeM771n_5-h(klA4R@Qg4lD?}I^*ybiskOYO(Q=wr%W66;qwi}e GP56IDoAucM diff --git a/examples/gringarten_multiapp/mrect1.i b/examples/gringarten_multiapp/mrect1.i index 6d3cee4a..15528d9e 100644 --- a/examples/gringarten_multiapp/mrect1.i +++ b/examples/gringarten_multiapp/mrect1.i @@ -1,11 +1,51 @@ -endTime = 3.16e8 - +endTime = 3e8 +frac_half_space = 20 #50 +n_elem_x = 7 #10 [Mesh] - uniform_refine = 0 - [matrixmesh] - type = FileMeshGenerator # reading mesh from file - file = 'gringartenbrickdomain100by100by10.e' #mesh created in cubit + # [matrixmesh] + # type = FileMeshGenerator # reading mesh from file + # file = 'gringartenbrickdomain100by100by10.e' #mesh created in cubit + # [] + [left_side] + type = GeneratedMeshGenerator + dim = 3 + nx = ${n_elem_x} + xmin = '${fparse -frac_half_space}' + xmax = 0 + bias_x = 0.8 + ny = 50 + ymin = -50 + ymax = 50 + nz = 1 + zmin = -5 + zmax = 5 + [] + [right_side] + type = GeneratedMeshGenerator + dim = 3 + nx = 7 + xmin = 0 + xmax = '${fparse frac_half_space}' + bias_x = 1.2 + ny = 50 + ymin = -50 + ymax = 50 + nz = 1 + zmin = -5 + zmax = 5 + [] + [right_side_rename] + type = RenameBoundaryGenerator + input = right_side + old_boundary = 'top bottom front back' + new_boundary = 'rtop rbottom rfront rback' + [] + [smg] + type = StitchedMeshGenerator + inputs = 'left_side right_side_rename' + clear_stitched_boundary_ids = true + stitch_boundaries_pairs = 'right left' [] [] @@ -39,11 +79,11 @@ endTime = 3.16e8 [Functions] [insitu_pp] type = ParsedFunction - value = '9.81*1000*(3000 - z)' + expression = '9.81*1000*(3000 - z)' [] [insitu_T] type = ParsedFunction - value = '363' + expression = '363' [] [] @@ -84,18 +124,26 @@ endTime = 3.16e8 [] [Preconditioning] - [./superlu] + [superlu] type = SMP full = true petsc_options_iname = '-ksp_type -pc_type -pc_factor_mat_solver_package' petsc_options_value = 'gmres lu superlu_dist' - [../] + [] [] +[Functions] + [dts] + type = ParsedFunction + expression = if(t<1e6,t*t,1e6) + [] +[] [Executioner] type = Transient solve_type = NEWTON - dt = 0.1e7 + # dt = 1e6 + dtmin = 1 + dtmax = 1e6 end_time = ${endTime} line_search = 'none' automatic_scaling = true @@ -105,6 +153,11 @@ endTime = 3.16e8 nl_max_its = 20 nl_rel_tol = 5e-04 nl_abs_tol = 1e-09 + [TimeStepper] + type = FunctionDT + function = dts + min_dt = 2 + [] [] [Outputs] @@ -184,50 +237,43 @@ endTime = 3.16e8 [Transfers] [normal_x_from_fracture] type = MultiAppNearestNodeTransfer - direction = from_multiapp - multi_app = fracture_app + from_multi_app = fracture_app source_variable = normal_dirn_x variable = fracture_normal_x [] [normal_y_from_fracture] type = MultiAppNearestNodeTransfer - direction = from_multiapp - multi_app = fracture_app + from_multi_app = fracture_app source_variable = normal_dirn_y variable = fracture_normal_y [] [normal_z_from_fracture] type = MultiAppNearestNodeTransfer - direction = from_multiapp - multi_app = fracture_app + from_multi_app = fracture_app source_variable = normal_dirn_z variable = fracture_normal_z [] [element_normal_length_to_fracture] type = MultiAppNearestNodeTransfer - direction = to_multiapp - multi_app = fracture_app + to_multi_app = fracture_app source_variable = element_normal_length variable = enclosing_element_normal_length [] [element_normal_thermal_cond_to_fracture] type = MultiAppNearestNodeTransfer - direction = to_multiapp - multi_app = fracture_app + to_multi_app = fracture_app source_variable = normal_thermal_conductivity variable = enclosing_element_normal_thermal_cond [] [T_to_fracture] - type = MultiAppInterpolationTransfer - direction = to_multiapp - multi_app = fracture_app + type = MultiAppGeometricInterpolationTransfer + to_multi_app = fracture_app source_variable = matrix_T variable = transferred_matrix_T [] [heat_from_fracture] type = MultiAppReporterTransfer - direction = from_multiapp - multi_app = fracture_app + from_multi_app = fracture_app from_reporters = 'heat_transfer_rate/joules_per_s heat_transfer_rate/x heat_transfer_rate/y heat_transfer_rate/z' to_reporters = 'heat_transfer_rate/transferred_joules_per_s heat_transfer_rate/x heat_transfer_rate/y heat_transfer_rate/z' [] diff --git a/examples/gringarten_multiapp/mrect1_amr.i b/examples/gringarten_multiapp/mrect1_amr.i deleted file mode 100644 index d0d015b3..00000000 --- a/examples/gringarten_multiapp/mrect1_amr.i +++ /dev/null @@ -1,265 +0,0 @@ -endTime = 3.16e8 - -[Mesh] - uniform_refine = 0 - [generate] - type = GeneratedMeshGenerator - dim = 3 - nx = 10 - xmin = -50 - xmax = 50 - ny = 10 - ymin = -50 - ymax = 50 - nz = 1 - zmin = -5 - zmax = 5 - [] -[] - -[GlobalParams] - PorousFlowDictator = dictator -[] - -[Variables] - [matrix_P] - [] - [matrix_T] - [] -[] - -[ICs] - [matrix_P] - type = FunctionIC - variable = matrix_P - function = insitu_pp - [] - [matrix_T] - type = FunctionIC - variable = matrix_T - function = insitu_T - [] -[] - -[BCs] -[] - -[Functions] - [insitu_pp] - type = ParsedFunction - value = '9.81*1000*(3000 - z)' - [] - [insitu_T] - type = ParsedFunction - value = '363' - [] -[] - -[PorousFlowFullySaturated] - coupling_type = ThermoHydro - porepressure = matrix_P - temperature = matrix_T - fp = water - gravity = '0 0 -9.81' - pressure_unit = Pa -[] - -[FluidProperties] - [water] - type = SimpleFluidProperties - thermal_expansion = 0 - [] -[] - -[Materials] - [porosity] - type = PorousFlowPorosityConst - porosity = 0.1 - [] - [permeability] - type = PorousFlowPermeabilityConst - permeability = '1E-16 0 0 0 1E-16 0 0 0 1E-16' - [] - [internal_energy] - type = PorousFlowMatrixInternalEnergy - density = 2875 - specific_heat_capacity = 825 - [] - [aq_thermal_conductivity] - type = PorousFlowThermalConductivityIdeal - dry_thermal_conductivity = '2.83 0 0 0 2.83 0 0 0 2.83' - [] -[] - -[Preconditioning] - [./superlu] - type = SMP - full = true - petsc_options_iname = '-ksp_type -pc_type -pc_factor_mat_solver_package' - petsc_options_value = 'gmres lu superlu_dist' - [../] -[] - -[Executioner] - type = Transient - solve_type = NEWTON - [TimeStepper] - type = IterationAdaptiveDT - dt = 100 - growth_factor = 2 - optimal_iterations = 20 - [] - dtmin = 1 - dtmax = 0.1e7 - end_time = ${endTime} - line_search = 'none' - # automatic_scaling = true - l_max_its = 20 - l_tol = 8e-3 - nl_forced_its = 1 - nl_max_its = 20 - nl_rel_tol = 5e-04 - nl_abs_tol = 1e-09 -[] - -[Adaptivity] - marker = points - max_h_level = 3 - stop_time = 1000 - [Markers] - [points] - type = ReporterPointMarker - x_coord_name = heat_transfer_rate/x - y_coord_name = heat_transfer_rate/y - z_coord_name = heat_transfer_rate/z - inside = refine - empty = do_nothing - [] - [] -[] - -[Outputs] - print_linear_residuals = false -[] - -[DiracKernels] - [heat_from_fracture] - type = ReporterPointSource - variable = matrix_T - value_name = heat_transfer_rate/transferred_joules_per_s - x_coord_name = heat_transfer_rate/x - y_coord_name = heat_transfer_rate/y - z_coord_name = heat_transfer_rate/z - [] -[] - -[Reporters] - [heat_transfer_rate] - type = ConstantReporter - real_vector_names = 'transferred_joules_per_s x y z' - real_vector_values = '10000; 0; 0; 0' - outputs = none - [] -[] - -[AuxVariables] - [normal_thermal_conductivity] - family = MONOMIAL - order = CONSTANT - [] - [fracture_normal_x] - family = MONOMIAL - order = CONSTANT - initial_condition = 0 - [] - [fracture_normal_y] - family = MONOMIAL - order = CONSTANT - initial_condition = 1 - [] - [fracture_normal_z] - family = MONOMIAL - order = CONSTANT - initial_condition = 0 - [] - [element_normal_length] - family = MONOMIAL - order = CONSTANT - [] - [fracDensity_AMR] - [] -[] - -[AuxKernels] - [normal_thermal_conductivity_auxk] - type = ConstantAux - variable = normal_thermal_conductivity - value = 5 - [] - [element_normal_length_auxk] - type = PorousFlowElementLength - variable = element_normal_length - direction = 'fracture_normal_x fracture_normal_y fracture_normal_z' - [] -[] - -[MultiApps] - [fracture_app] - type = TransientMultiApp - input_files = frect1.i - execute_on = TIMESTEP_BEGIN - sub_cycling = true - [] -[] - -[Transfers] - [normal_x_from_fracture] - type = MultiAppNearestNodeTransfer - direction = from_multiapp - multi_app = fracture_app - source_variable = normal_dirn_x - variable = fracture_normal_x - [] - [normal_y_from_fracture] - type = MultiAppNearestNodeTransfer - direction = from_multiapp - multi_app = fracture_app - source_variable = normal_dirn_y - variable = fracture_normal_y - [] - [normal_z_from_fracture] - type = MultiAppNearestNodeTransfer - direction = from_multiapp - multi_app = fracture_app - source_variable = normal_dirn_z - variable = fracture_normal_z - [] - [element_normal_length_to_fracture] - type = MultiAppNearestNodeTransfer - direction = to_multiapp - multi_app = fracture_app - source_variable = element_normal_length - variable = enclosing_element_normal_length - [] - [element_normal_thermal_cond_to_fracture] - type = MultiAppNearestNodeTransfer - direction = to_multiapp - multi_app = fracture_app - source_variable = normal_thermal_conductivity - variable = enclosing_element_normal_thermal_cond - [] - [T_to_fracture] - type = MultiAppInterpolationTransfer - direction = to_multiapp - multi_app = fracture_app - source_variable = matrix_T - variable = transferred_matrix_T - [] - [heat_from_fracture] - type = MultiAppReporterTransfer - direction = from_multiapp - multi_app = fracture_app - from_reporters = 'heat_transfer_rate/joules_per_s heat_transfer_rate/x heat_transfer_rate/y heat_transfer_rate/z' - to_reporters = 'heat_transfer_rate/transferred_joules_per_s heat_transfer_rate/x heat_transfer_rate/y heat_transfer_rate/z' - [] -[] diff --git a/examples/gringarten_multiapp/plot_compare.py b/examples/gringarten_multiapp/plot_compare.py index 72f9c8ed..923e5fcc 100755 --- a/examples/gringarten_multiapp/plot_compare.py +++ b/examples/gringarten_multiapp/plot_compare.py @@ -1,31 +1,34 @@ #!/usr/bin/env python3 -#* This file is part of the MOOSE framework -#* https://www.mooseframework.org -#* -#* All rights reserved, see COPYRIGHT for full restrictions -#* https://github.com/idaholab/moose/blob/master/COPYRIGHT -#* -#* Licensed under LGPL 2.1, please see LICENSE for details -#* https://www.gnu.org/licenses/lgpl-2.1.html +# * This file is part of the MOOSE framework +# * https://www.mooseframework.org +# * +# * All rights reserved, see COPYRIGHT for full restrictions +# * https://github.com/idaholab/moose/blob/master/COPYRIGHT +# * +# * Licensed under LGPL 2.1, please see LICENSE for details +# * https://www.gnu.org/licenses/lgpl-2.1.html import pandas as pd -import os -import sys import numpy as np import matplotlib.pyplot as plt import math -#------------------------------------------------------------------------------- -df = pd.read_csv('mrect1_out_fracture_app0.csv') -print('fracture dataframe description: ',df.columns.values) -['time' 'J_out' 'P_in' 'P_out' 'TK_in' 'TK_out' 'kg_out' 'kg_per_s'] -df_amr = pd.read_csv('mrect1_amr_out_fracture_app0.csv') -#analytic solution from Koenradd and Gringarten paper, 1975 -timevector = np.logspace(1, 9, 50, base=10) -print(timevector) - -#analytic Gringarten solution 1975 from Beckers, Koenraad -Trock = 363-273.15 -Tinj = 303-273.15 +from mpmath import * + + +# ------------------------------------------------------------------------------- +df = pd.read_csv("mrect1_out_fracture_app0.csv") +print("fracture dataframe description: ", df.columns.values) +["time" "J_out" "P_in" "P_out" "TK_in" "TK_out" "kg_out" "kg_per_s"] +df_50 = pd.read_csv("mrect1_out_fracture_app0_d50.csv") +# analytic solution from Koenradd and Gringarten paper, 1975 +timevector = np.logspace(6, 11, 50, base=10) + +# analytic solution for single fracture Gringarten 1975 eqn A18 +# from Beckers, Koenraad + +# terms for infinite and finite space fractures +Trock = 363 - 273.15 +Tinj = 303 - 273.15 rhowater = 1000 cpwater = 4150 fracnumb = 1 @@ -35,21 +38,86 @@ fracheight = 10 rhorock = 2875 cprock = 825 -mtot = 25/25/10 -Q = mtot/rhowater/fracnumb/fracwidth -td = (rhowater*cpwater)*(rhowater*cpwater)/(4*krock*rhorock*cprock)*(Q/fracheight)*(Q/fracheight)*timevector -Twater=[] -for i in range(len(timevector)): - temp=1-math.erf(1/(2*math.sqrt(td[i]))) - Twater.append(Trock-temp*(Trock-Tinj)) - - -#------------------------------------------------------------------------------- -fig1=plt.figure() -plt.plot(timevector/3600/24/365,Twater,'b', linestyle = '-',marker = 'o',fillstyle='none',label='Gringarten 1975') -plt.plot(df['time']/3600/24/365,df['TK_out']-273.15,'r', linestyle = '--',label='Falcon') -plt.plot(df_amr['time']/3600/24/365,df_amr['TK_out']-273.15,'g', linestyle = '-.',label='Falcon amr') -plt.ylim([50, 100]) +mtot = 25 / 25 / 10 +Q = mtot / rhowater / fracnumb / fracwidth +# terms needed for finite spaced fractures. Gringarten 1975 eqn A19 +xe = 20 # 1/2 spacing between fractures +z = fracheight + + +# nonDimensionalized Time eqn 10 +td = ( + (rhowater * cpwater) + * (rhowater * cpwater) + / (4 * krock * rhorock * cprock) # where does the 4 come from? Its not in eqn 10 + * (Q / fracheight) + * (Q / fracheight) + * timevector +) +# Solution to eqn A18 -- infinite space between fractures +Twater = [] +Twd_inf = [] +for i in range(len(td)): + # eqn A18 solving for nonDimensionalized temperature + temp = 1 - math.erf(1 / (2 * math.sqrt(td[i]))) + Twd_inf.append(temp) + # convert to Temperature + Twater.append(Trock - temp * (Trock - Tinj)) + +# Solution to eqn A19 -- finite space between fractures set above by xe +# simplified inverse laplace transform +mp.dps = 15 +top = rhowater * cpwater * Q * xe +bottom = 2 * krock * fracheight +zd = z / fracheight +Twd_tilde = lambda s: 1 / s * exp(-zd * sqrt(s) * tanh(top / bottom * sqrt(s))) +Twater_finite = [] +for i in range(len(td)): + temp = invertlaplace(Twd_tilde, td[i], method="talbot") + # convert to Temperature + Twater_finite.append(Trock - temp * (Trock - Tinj)) + +# ------------------------------------------------------------------------------- +fig1 = plt.figure() +plt.plot( + timevector / 3600 / 24 / 365, + Twater, + "b", + linestyle="-", + linewidth=1, + marker="o", + alpha=0.5, + fillstyle="none", + label="Gringarten xe=inf", +) +plt.plot( + timevector / 3600 / 24 / 365, + Twater_finite, + "k", + linestyle="-", + linewidth=1, + marker="^", + alpha=0.5, + fillstyle="none", + label="Gringarten xe=20m", +) +plt.plot( + df_50["time"] / 3600 / 24 / 365, + df_50["TK_out"] - 273.15, + "r", + linestyle="--", + linewidth=2, + label="Falcon xe=50m", +) +plt.plot( + df["time"] / 3600 / 24 / 365, + df["TK_out"] - 273.15, + "g", + linestyle="-.", + linewidth=2, + label="Falcon xe=20m", +) +plt.ylim([30, 100]) plt.xlim([0, 10]) plt.ylabel("T (degC)") plt.xlabel("Time (year)") diff --git a/examples/gringarten_multiapp/production_single_fracture_amr.xyz b/examples/gringarten_multiapp/production_single_fracture_amr.xyz deleted file mode 100644 index 83008e41..00000000 --- a/examples/gringarten_multiapp/production_single_fracture_amr.xyz +++ /dev/null @@ -1 +0,0 @@ -0.01 0.1 50 0 diff --git a/examples/gringarten_multiapp/rect_100_10_medium.e b/examples/gringarten_multiapp/rect_100_10_medium.e deleted file mode 100644 index 2a3ae8c75d453911f99c521e6788ecbc0803629b..0000000000000000000000000000000000000000 GIT binary patch literal 0 HcmV?d00001 literal 14484 zcmeI&hkH-=-^cN9L_*A<_6jOSh!sJE9E6lst9F$pK@d$OBKAmZT1BnWqA^0$ri%V( zP8GGQt+rYnX^ZaCs?njjpRbd1PQJbFtJnPp{O<2{J@WpH_ZgpYjtpI*yA|;id`b{X zqIgomnAkD#$q62hk9q>q6O$8S$E7Evs#fqRHfDTs?C`{7)qD}4`k0jA3F8X-&l8HP zEg>m^Bp-u%iGq6U(4?oWB%QpB-Ih2;%`Mb!<(9S27n_=p#yPAGUB^1d*!bAAgkj3H zZm{Txl<{eli=u??eO!84;+T=@_S8Ws6hq1Q)Wp~c32EaJQ^rvB_`#=>XUdnvV*k=! zWkga+d^%@~O-)H0lO8)UC2e@@xWp;UjdeYZp*R9YBqlLO1sz^}>aUO-KPfR;^TV9R za)Tq&N73vaA4ItnCl3$f&g0X^YYH_M9jlLDdSZH#F33k~)n(ZDp^51YyF@e_5uY?H zWlXbVW+S~tGkuM1(Y!fglM{v~j!zCvXcQ6Dyj`;vVa-}LZ`mRus71?&u<(eM;numd z)}>LxTHlGIQtXT9y3Vzg=nsvNPumrI2}l@Pc%#&ntGRSOx_41n%a`Khr>3Q(Y7w1W zt%=G-3-*~cF5m9av=kic9G^wqe^?5$K1p5HG1>>thkNt4CX>T9(^!x@+3ME$|0TDi zu3wNl<^RepmrjPppj*}2Xbrr_s}s7sy7#s?h9PQTE<%$HyOuj8_3w3@T}%Rz2?_KeWMUGJATl0Bp8 zq2sHsd$l)RdUS1V@A|^|v*vi$7ah3qt#|!S(S;tbdhp)xLwnr1bKi`&NUtgxZade1 zrFzPCuX>o@fmTm?)k8O5{dJvJJ$&9ZSM%w3bqqEK{PrlbPZ(Fdt zPn^7*ufF|#V(*V5bNq(is2trrGN)7U+*S3DN9J@mIpLX;m)!H@bm~8K;pnnn^#Qv# z>?pimPN%-1dul$(@qNw(c3t3|clz>`rnC?EzSRTOz1sWxUV8aM+WU4qTIr5g`+%`k zJALKUwJrm0^|*Azz28p#D~#wkRL4b>?0VKKujzrRw-$Tl4fuM?H*3B7SHvv$K2N_1 zqw19AJpId7?Q_G)quyt2c=S=l{sTq_Iqlp>(1$41_Yv^=$)!#{ z^@ja|S3Or)pZ@3gPdD%Ns{7Ab`B08iXI-BOIDhxu$Q6C3D`f60IslrXAW`Gv)9w;>I}*4!=z7HOup5En?*2H4%f{{Io;u z+GUT~M0=Y?%U5h6)_?Ls>my$MZ9aat*IL?Zzc_bZQ{wc!^DpNTyNw@t>^rZ#sqK63 zt3rFLdlioUT-bi_?wJWg0t@Sh7QMHla#!+Nr*-OopBOPOw2psa|H11+&eTh%{iUjp z?s`Y(Tk%DoQ^bSEw$1W7@8M%L&BNCs)0X`db86K6$kYYb7POmmDKhQg-JhpyKgYTZ z+}-c9!R~q8`Rdev!q4BFblX|)qIW;|#I2_t%o}~X)1TTGzrD4;a2}}NhpSE6!TN4{PN(iQM@`6a>&*Y%c9X*yPrL7X?95cmD1B^IWLoOo zS9UL-uX^9;(7D96?SEV|gxG9r@AjpL{#*X&@h0P|`1>a(J~#Q-M<4S6JC@GBlDfK! zVRSF38!=y>eYk#VbYcDT;yeEsJGiiJT8aN&5124$`RP(d*U#psuN*MR=$2WKA8G}U z(cJLed9dGuMhm+=B(rtIFu;7p_UXUj8t-Hk717^>NlU`?W#q2410XJxANn*pW0v+WHVm-JZ7DY8FWSNDhSS!A}i+dJ9O_~OsV&)ofrGchB-)5EkA5S+9$@^dDbFLHf?S&dX zaedu>=RRVhX2w7Me7M>5n}&Z4KND$oiL1S}_OO9w^o9zRKiGcN9Z%<=BFB0^m+Fq! z_#*$Kzcvjvoq3AP%?~{MRkRrumH%LGpGWNb_hG*Cs|6nUk$KK9AM$LTyACaL^V>Cy z%VB)^^(HsE>5i{?RDS(IA60SI(`cFdQ-wuue``4756<{|!ln4tjQ{T*pUg0v@f)q9 zzFEEf)Fi`M*E$)051ZFI$k{hWo80SLyB$b1oPAI_|Bqe^!aU5|{kBzZ-*ld*sNXzg zuhdX^kh!|o9@03&MXprBqYy$E~28Y&rV;k$FTR85t8Y1x%S35_c>~AcsAXf zlxH~4pYg@5$0z)C(Xi*;=v@2OyQ9X38TLLj*q7s@haNHPeP~3#eB3|p+(5(LhemW( ztBm(YryE`3is#OIY=@f{Q*qW1^UkM+voDR@>)R z%i=xoz_qt-earnuG!MVAhjB%QM{Zo|JlBkS@rto|K$$hfAn&n;}#j$PCfI+J+o>4g(6|y@43%! ztNd%b?)jZ_pWC+i*B)H9)nkTdT>tLgN0-g^8DA%zu9f9Jza8?geo{8@h}kLQ*DK?? zuQ!dTTd^76>@qv&{yzQu@@L%VHQJN9-xF}v>>7RfcG7@*X1DdTMlag&qS-M2Y--{S z(`adW29LZq(0vcu`_pJ`c!EFL7U%3kqmAjQe74*vXWvD&?>273H0>wDQ?X>TICnoq znx2XcKU!e#ACD1bc&cSIiF5Zwn5W#u*>$oEBhS;|@!dw@`Um+`998#Kcm2aN1Mj?_ z6v*`gFMjjsB+egrVL;ABqife{v88^v%J~DY-Z*fR^A^9n?m`vLQ!4oBI=wi)e#-fn zw(N%w#?~tl=C0oZ1@9rlcUq0Pi)R`IHIFCv?!!M6E17CI>lk%=cGQ8^@I;{9deeHohL0)Gwryr%z31C@$br0?JxOvp0-+j(KF6tN0 z(9auh_i+8dE@y@r-g)oD?8`Lm^SIGIf9|cLza)jfU^@H6?e|}@W#DHQf=s`Ymri7i z%`^S(`;^FUX_|h&|8j8Jio)+_{O&iIdwf}*dz|0>rt2z9F1()K{l+Ob{j)f}ai#e! z4DY;GKYB5f_NMQJ{ONssld4grZH0U9D-RSAWQwfF8QvCiNz?%X9LmH~OLen|X@cFVz3X zeHekO^1P4V*>RJ5oUY$7JJ!3uN%gYptn>Sf*_Y*Y{w7Q7H!z*LjteV&!uz<;bxRJ# zIsM#6^GEz!1UYr`E1t}V4KfP%OSxGSw`CR9r+#~B>Q0YWy~N3)B{Q8m>-zf5GIN-} zzs|?Ryy)JC>Hqzfu5lOgOlSU_==au~p0#Fqoc7#vOZ3yhoRZvaPhtMWdlN?$J)!pT3%h@j_3m%dqx>edU#y!o z!a84%NGpC{E2*Ma-tSjd%}${I!>vpAto7uc|Kr@&C#EPrGJJishjj^^^Yyq)XI|XM zeww+uX<<9{YB8nroPLfEzBX^Ej&D+ORmMN^{9kQtI(hVmTwnO5;nX>9QJZNqbzI{^ zTmQ~-+NlSPjK3S?)XAH2E+|&_IjMNVcN8aXx%o5k;OFhTHnHlN!DkdF%<|+aj(L4k zxMKGuMRpPIraZN^;1`?DyW+nw7W~Tizj5RRW1jrqJ9KO0|NjJ7zdtX+WaxM7`YuAh zyw~4zirHngdH+p7e-pC)HbM-5{ywDKGIm*iU38$qapx~ONDgNP5JKJ@$C26kC(4WXaZ8zaOntF4J$Z%Wip zc+JopE$p(|TH5ti#4x+wny9~!son-{?Rq<61ll7K9iX2`bbb;{59)|iJvCC@fiAV7mdfD~IiM{Q5AENF@_o8w3ML+b%01U(v z(0%G2HD819BnD%MU4M!gi#TX~<1rM&FdPYp!w9HvBs4CKNozI=iBR5XJdGqIV+>M| zim^~$8q_xq=@<|7O~6Dv5inX;9vDsBZ=`?0P0~CT2l-voXi6&n3>Y>+^{V z?D|6DGk6vyu?WvW_r4gq{{R$22`s@DiQ zC|Y3$)Tg$VP+My}ik)~2yP!Pvhe7>q5N^kI#P(>5NbH7=QLft5ul7B77>}R>bYHsX zC}`XUl%w_gO?^@iqDQEi%c)it(e&|0WZ?OKbz(0V9peYD0}Yt>aA z0Ij+9gw|a1(g@lo>QnnB3__M2wRfJv5IbuBX|HLH>SHL>r*_R%JTzaiD35a3hy8dJ z>K}%YP~He=uMbBYbld?b_cf?rxr?E^r}2^<)u(pVM`JBUApttyLFgQZp#CK2c~f33 zXf4M;eQH098c<&{QZNo{kcQPrMIv;a*P(O&2kKYuDx_m9lzRl~KZ=P^o^mE&5~g4! zCS!#i)i%SfPbV(JRE&qt{RVWtT&Q2UnNVH^mfJCtsCid?A56@|H0XS9LgzRJ_3K#W z&&Dh$_bsUZI810>W1%(DFQ`A&^ZYinUdmHhQT-|_YK^q+T5IL0tf+pK6}9dspuMC#l@-;mvf?P1&^}b2%8Ke& zSy6jZ`(JyuoL$yF*D-20p}u9%v!p(im0JNTpl4HkDy!`Qtc0Fn^{K43N_YVip+1!r z)vvPpt6(+Kpgffo)vvPpt3mIaWGGK%wbj63jDh-8R#d;ripp7w(NLbsit1Ncu@+v! z(@>tuzKF(pJJzwwBZ$gVyYiH?4)IW*%F3;Wy%-AhsjRln*nlBWpUR5rS6NXxSr`Q6 zsjR4emFuGcHsMJqPi3`5VG9O8eJU%eUu8w*WMd$dr?R5@Rc?e3Y(-xvPh~~*tE~Pe z*oJ;kp2}+LfbHlF^{K3=ew7uK^9p)Fc`7$U3+zBos83}@^{cG@me_?yp*)op)vvPp z!>}8VKzS;wEgXB$7M<`AV(>83u6*UHtQg*qYTO-5am!F&huK4vYx-nsDi4{GgTdWW@QE!se9(GC%4k4R{LbcFV* z0j-sui!RVJ*A?BMwbL4QM-ONnwT9XcTKC7$3y-5W`asL6{n#Jc_uAJ_K+ln$mBAQ- zrx1%c#A7IiVK@>n0wXaBi5QKik%VN7K?+ha7HJrVbd1LYOvEHi#uQA&G)%_~WFQkW zF$=RX2XiqG^RWO6@eH2DB0Pu1Sc0W^9?P&CE3gu)@B&t24PL}cSc`R7kC%~!4cLfH z*o@GklIOkcThv6~4w*T*GyI zgKzO2zQ+&v5kKK)+`vuzf?sh9zu`9S;4bdrcl?1r@fYsnZ{*wmebt9p1imN=Klq~< ziX#9eP!gq38f8!xfhdRasDKAh5tUFGRZtbxP#rZ;6SYtqL8ybesD}qpAHisVhG>Mw z2tgAxMJSq~Ia;74S|JRr5so%!i*|@WdqkoGIwA^&hc96K|Ns7D4x-T&-4KI^&>cPS zFdjiqJc`HA3y-5W`k*iRp+5#-AfCV=Jc+>=f~OFRIK*QphG94oFajen3W*qvr;&tY zj6n)gF&1eUhjfg`1Wd#vOvV&U#WYOE3}he^GcgOZF$Z%o5A(4A3-Jt|#Ueb1#aM!+ zcpl5J94oLAtMCF=V+~%!OIV9_SdW*Hg$>wjb+cI?1T?80vB!Cvgc ze!PkUcnt?}2#4`H{(~bpiZ_r86K~=e-okO5z)8G~Q#g$?IE!;Qj|+GQ7x6CM!~3{| z5AY#A!exAnPjCgF;xl}XFOY{X@fE(tRb0b$e1mWC9lpm8_z^$hXWYO|{DNO`3%}tu z?%*!&;dlIjKk*mt<8R~_Jo-L9#3JxTQTV|h#ZVjpD1nkFh0-X4vIs;ult%?TfQqPu z%BX^>sD|pOftsj=+6Y1&)I~i!i24Xd12jYR?-%CZdcUw|+~4yXlox`% zWW8@X5f_2|WxaErC3e9QEJjx>1$)qX2X!MZ1N+i?4=pD?gq7%yRbby*@1-8Z)p!_d zz@E0=Pwa2&-Bgme0qlGCXSaTob1(p%u@(KX4F~WlUc*7e;1CYubwuMRjvyCrfc94DRx``3C8ohF{b8L-Fy@h&o{zlcZB6KnA(*5NVy2QQ-+vhX-I cqBr(q6Z&8?`eGmUVh>)yZtTKN?7;T_0RX(KtN;K2 diff --git a/examples/gringarten_multiapp/tests b/examples/gringarten_multiapp/tests index da863bd3..3a6aa910 100644 --- a/examples/gringarten_multiapp/tests +++ b/examples/gringarten_multiapp/tests @@ -1,11 +1,10 @@ [Tests] [shortTest] - method = '!dbg' type = CSVDiff input = 'mrect1.i' csvdiff = 'mrect1_out_fracture_app0.csv' - cli_args = 'endTime=6e7' - requirement = 'Run Gringarten multiapp for a few timesteps.' - skip = 'long time' + cli_args = 'endTime=6e6' +requirement = 'The system shall match the Gringarten analytic solution for a set of parallel fracturs seperated by a finite distance using a multi-app to transfer heat between a 2D fracture and 3D matrix.' + heavy = true [../] []