.. DO NOT EDIT. .. THIS FILE WAS AUTOMATICALLY GENERATED BY SPHINX-GALLERY. .. TO MAKE CHANGES, EDIT THE SOURCE PYTHON FILE: .. "auto_examples_mesh/example_meshconvergence_thacker3d.py" .. LINE NUMBERS ARE GIVEN BELOW. .. only:: html .. note:: :class: sphx-glr-download-link-note Click :ref:`here ` to download the full example code .. rst-class:: sphx-glr-example-title .. _sphx_glr_auto_examples_mesh_example_meshconvergence_thacker3d.py: ================================================================================ Thacker3d mesh convergence ================================================================================ In this example mesh convergence is carried out in the case of the thacker3d with variable density. .. GENERATED FROM PYTHON SOURCE LINES 10-25 .. code-block:: default import os import h5py import time import numpy as np import matplotlib.pyplot as plt from matplotlib import cm import matplotlib.image as mpimg from mpl_toolkits.axes_grid1 import make_axes_locatable, axes_size import freshkiss3d as fk import freshkiss3d.extra.plots as fk_plt #sphinx_gallery_thumbnail_number = 6 os.system('rm -r outputs') .. rst-class:: sphx-glr-script-out .. code-block:: none 256 .. GENERATED FROM PYTHON SOURCE LINES 26-28 Mesh parameters: -------------------- .. GENERATED FROM PYTHON SOURCE LINES 29-35 .. code-block:: default N_mesh = 5 #Number of mesh for mesh convergence First_mesh = 2 NL_start = 10 NL = [NL_start+2*r for r in range(N_mesh)] .. GENERATED FROM PYTHON SOURCE LINES 36-38 Case parameters: -------------------- .. GENERATED FROM PYTHON SOURCE LINES 39-63 .. code-block:: default H0 = 0.1 ETA = 0.1 ALPHA = 1. A = 10. L = 4. OMEGA = np.sqrt(ALPHA*9.81) PERIOD = 2.*np.pi/OMEGA RHO_0 = 998. FINAL_TIME = 1.*PERIOD #Mesh options: MESH_SPLITTING = False #Plot options: PLOT_MESH = True PLOT_ERROR = True #Error options: DATA_REP = "data" ERROR_TYPE = 'L2' SAVE_ERROR = True SECOND_ORDER = True .. GENERATED FROM PYTHON SOURCE LINES 64-66 Lists initialization: ----------------------- .. GENERATED FROM PYTHON SOURCE LINES 67-78 .. code-block:: default cases = [] simutime_list = [] triangular_mesh_list = [] NT = [] layer_list = [] primitives_list = [] problem_list = [] thacker3d_analytic_list = [] CPUTime_list = [] .. GENERATED FROM PYTHON SOURCE LINES 79-81 Meshes: -------------------- .. GENERATED FROM PYTHON SOURCE LINES 82-118 .. code-block:: default os.system('gmsh -2 -format msh2 inputs/thacker2d{}.geo -o inputs/thacker2d{}.msh'.format(First_mesh, First_mesh)) TG, vertex_labels, boundary_labels = fk.read_msh('inputs/thacker2d{}.msh'.format(First_mesh)) os.system('rm inputs/thacker2d{}.msh'.format(First_mesh)) x = np.asarray(TG.x) y = np.asarray(TG.y) trivtx = np.asarray(TG.trivtx) x *= 0.12 y *= 0.12 x -= 0.6 y -= 0.6 triangular_mesh_list.append( fk.TriangularMesh(TG, vertex_labels, boundary_labels) ) NT.append(triangular_mesh_list[0].NT) if PLOT_MESH: fk_plt.plot_mesh(triangular_mesh_list[0]) for M in range(1, N_mesh): if MESH_SPLITTING: triangular_mesh_list.append( triangular_mesh_list[M-1].refine_by_splitting() ) NT.append( triangular_mesh_list[M].NT) if PLOT_MESH: fk_plt.plot_mesh(triangular_mesh_list[M]) else: mesh_id = First_mesh+M os.system('gmsh -2 -format msh2 inputs/thacker2d{}.geo -o inputs/thacker2d{}.msh'.format(mesh_id, mesh_id)) TG, vertex_labels, boundary_labels = fk.read_msh('inputs/thacker2d{}.msh'.format(mesh_id)) os.system('rm inputs/thacker2d{}.msh'.format(mesh_id)) x = np.asarray(TG.x) y = np.asarray(TG.y) trivtx = np.asarray(TG.trivtx) x *= 0.12 y *= 0.12 x -= 0.6 y -= 0.6 triangular_mesh_list.append( fk.TriangularMesh(TG, vertex_labels, boundary_labels) ) NT.append(triangular_mesh_list[M].NT) if PLOT_MESH: fk_plt.plot_mesh(triangular_mesh_list[M]) .. rst-class:: sphx-glr-horizontal * .. image-sg:: /auto_examples_mesh/images/sphx_glr_example_meshconvergence_thacker3d_001.png :alt: Mesh :srcset: /auto_examples_mesh/images/sphx_glr_example_meshconvergence_thacker3d_001.png :class: sphx-glr-multi-img * .. image-sg:: /auto_examples_mesh/images/sphx_glr_example_meshconvergence_thacker3d_002.png :alt: Mesh :srcset: /auto_examples_mesh/images/sphx_glr_example_meshconvergence_thacker3d_002.png :class: sphx-glr-multi-img * .. image-sg:: /auto_examples_mesh/images/sphx_glr_example_meshconvergence_thacker3d_003.png :alt: Mesh :srcset: /auto_examples_mesh/images/sphx_glr_example_meshconvergence_thacker3d_003.png :class: sphx-glr-multi-img * .. image-sg:: /auto_examples_mesh/images/sphx_glr_example_meshconvergence_thacker3d_004.png :alt: Mesh :srcset: /auto_examples_mesh/images/sphx_glr_example_meshconvergence_thacker3d_004.png :class: sphx-glr-multi-img * .. image-sg:: /auto_examples_mesh/images/sphx_glr_example_meshconvergence_thacker3d_005.png :alt: Mesh :srcset: /auto_examples_mesh/images/sphx_glr_example_meshconvergence_thacker3d_005.png :class: sphx-glr-multi-img .. GENERATED FROM PYTHON SOURCE LINES 119-124 Cases: -------------------- Here we define cases based on a ``case`` class containing info such as space order or time order for each simulation we want to run. .. GENERATED FROM PYTHON SOURCE LINES 125-143 .. code-block:: default class case(): def __init__(self, NL=1, mesh_id=0, time_order=False, space_order=False): self.NL = NL self.mesh_id = mesh_id self.time_order = time_order self.space_order = space_order self.ipres = True if space_order is True or time_order is True: self.order = 2 else: self.order = 1 for M in range(N_mesh): cases.append(case(NL=NL[M], mesh_id=M, time_order=False, space_order=False)) if SECOND_ORDER: cases.append(case(NL=NL[M], mesh_id=M, time_order=True, space_order=True)) .. GENERATED FROM PYTHON SOURCE LINES 144-146 Boundary conditions -------------------- .. GENERATED FROM PYTHON SOURCE LINES 147-150 .. code-block:: default fluvial_heights = [ fk.FluvialHeight(ref=r, height=0.0) for r in [1,2,3,4] ] .. GENERATED FROM PYTHON SOURCE LINES 151-153 Solving (loop over all cases): ------------------------------- .. GENERATED FROM PYTHON SOURCE LINES 154-261 .. code-block:: default for I, case in enumerate(cases): start = time.time() print(" ") print(" ") print(" SOLVING CASE {} ".format(I)) print(" ") print(" ") mesh_id = case.mesh_id time_second_order = case.time_order NUM_PARAMS={'ipres':case.ipres, 'space_second_order':case.space_order, 'velocity_correction':False, 'rho_correction':True, 'flux_type':2} PHY_PARAMS={'variable_density':True, 'variable_density_model':'ns-fourier'} #'ns-fourier' or 'boussinesq' simutime_list.append(fk.SimuTime( final_time=FINAL_TIME, time_iteration_max=10000, second_order=time_second_order)) if PHY_PARAMS['variable_density_model'] == 'boussinesq': HAS_TRACER = True HAS_DENSITY = False elif PHY_PARAMS['variable_density_model'] == 'ns-fourier': HAS_TRACER = False HAS_DENSITY = True thacker3d_analytic_list.append(fk.Thacker3DVariableDensity( triangular_mesh_list[mesh_id], NL=case.NL, a=ALPHA, h0=H0, eta=ETA, L=L, b=A, density = HAS_DENSITY, tracer = HAS_TRACER, compute_error=True, error_type=ERROR_TYPE, error_output='none')) thacker3d_analytic_list[I](0.) layer_list.append(fk.Layer( case.NL, triangular_mesh_list[mesh_id], topography=thacker3d_analytic_list[I].topography)) primitives_list.append(fk.Primitive( triangular_mesh_list[mesh_id], layer_list[I], height=thacker3d_analytic_list[I].H, Uinit=thacker3d_analytic_list[I].U, Vinit=thacker3d_analytic_list[I].V, HRhoinit=thacker3d_analytic_list[I].HRho)) problem_list.append(fk.Problem( simutime_list[I], triangular_mesh_list[mesh_id], layer_list[I], primitives_list[I], analytic_sol=thacker3d_analytic_list[I], numerical_parameters=NUM_PARAMS, physical_parameters=PHY_PARAMS, fluvial_heights=fluvial_heights)) problem_list[I].solve() CPUTime_list.append(time.time()-start) error_H = np.zeros((N_mesh, 2)) error_QX = np.zeros((N_mesh, 2)) error_QY = np.zeros((N_mesh, 2)) error_Rho = np.zeros((N_mesh, 2)) error_H_cumul = np.zeros((N_mesh, 2)) error_QX_cumul = np.zeros((N_mesh, 2)) error_QY_cumul = np.zeros((N_mesh, 2)) error_Rho_cumul = np.zeros((N_mesh, 2)) CPUTime = np.zeros((N_mesh, 2)) for I, case in enumerate(cases): CPUTime[case.mesh_id, case.order-1] = CPUTime_list[I] sol = problem_list[I].analytic_sol error_H[case.mesh_id, case.order-1] = sol.error_H error_QX[case.mesh_id, case.order-1] = sol.error_QX error_QY[case.mesh_id, case.order-1] = sol.error_QY error_Rho[case.mesh_id, case.order-1] = sol.error_HRho error_H_cumul[case.mesh_id, case.order-1] = sol.error_H_cumul/FINAL_TIME error_QX_cumul[case.mesh_id, case.order-1] = sol.error_QX_cumul/FINAL_TIME error_QY_cumul[case.mesh_id, case.order-1] = sol.error_QY_cumul/FINAL_TIME error_Rho_cumul[case.mesh_id, case.order-1] = sol.error_HRho_cumul/FINAL_TIME if SAVE_ERROR: os.system('mkdir {}'.format(DATA_REP)) os.system('rm {}/thacker3d_errors.h5'.format(DATA_REP)) REF_FILE = DATA_REP+"/thacker3d_errors.h5" with h5py.File(REF_FILE, 'w') as f: f.create_dataset('CPUTime', data=CPUTime) f.create_dataset('NT', data=NT) f.create_dataset('NL', data=NL) f.create_dataset('error_H', data=error_H) f.create_dataset('error_QX', data=error_QX) f.create_dataset('error_QY', data=error_QY) f.create_dataset('error_Rho', data=error_Rho) f.create_dataset('error_H_cumul', data=error_H_cumul) f.create_dataset('error_QX_cumul', data=error_QX_cumul) f.create_dataset('error_QY_cumul', data=error_QY_cumul) f.create_dataset('error_Rho_cumul', data=error_Rho_cumul) .. rst-class:: sphx-glr-script-out .. code-block:: none SOLVING CASE 0 =================================================================== | INITIALIZATION | =================================================================== Problem size: Nodes=463, Layers=10, Triangles=848, Iter = 0 ; Dt = 0.0000s ; Time = 0.00s ; ETA = 0.00s =================================================================== | TIME LOOP | =================================================================== Iter = 5 ; Dt = 0.0093s ; Time = 0.05s ; ETA = 2.23s Iter = 9 ; Dt = 0.0093s ; Time = 0.08s ; ETA = 1.32s Iter = 14 ; Dt = 0.0093s ; Time = 0.13s ; ETA = 1.22s Iter = 18 ; Dt = 0.0094s ; Time = 0.17s ; ETA = 1.21s Iter = 22 ; Dt = 0.0094s ; Time = 0.21s ; ETA = 1.14s Iter = 27 ; Dt = 0.0094s ; Time = 0.25s ; ETA = 1.14s Iter = 31 ; Dt = 0.0095s ; Time = 0.29s ; ETA = 1.09s Iter = 35 ; Dt = 0.0095s ; Time = 0.33s ; ETA = 1.08s Iter = 40 ; Dt = 0.0095s ; Time = 0.38s ; ETA = 1.00s Iter = 44 ; Dt = 0.0095s ; Time = 0.41s ; ETA = 1.02s Iter = 48 ; Dt = 0.0096s ; Time = 0.45s ; ETA = 0.98s Iter = 53 ; Dt = 0.0096s ; Time = 0.50s ; ETA = 0.98s Iter = 57 ; Dt = 0.0096s ; Time = 0.54s ; ETA = 0.93s Iter = 61 ; Dt = 0.0096s ; Time = 0.58s ; ETA = 0.93s Iter = 65 ; Dt = 0.0097s ; Time = 0.62s ; ETA = 0.88s Iter = 70 ; Dt = 0.0097s ; Time = 0.66s ; ETA = 0.86s Iter = 74 ; Dt = 0.0097s ; Time = 0.70s ; ETA = 0.82s Iter = 78 ; Dt = 0.0098s ; Time = 0.74s ; ETA = 0.83s Iter = 82 ; Dt = 0.0098s ; Time = 0.78s ; ETA = 0.78s Iter = 86 ; Dt = 0.0098s ; Time = 0.82s ; ETA = 0.77s Iter = 91 ; Dt = 0.0099s ; Time = 0.87s ; ETA = 0.73s Iter = 95 ; Dt = 0.0099s ; Time = 0.91s ; ETA = 0.71s Iter = 99 ; Dt = 0.0099s ; Time = 0.95s ; ETA = 0.66s Iter = 103 ; Dt = 0.0099s ; Time = 0.99s ; ETA = 0.66s Iter = 107 ; Dt = 0.0099s ; Time = 1.03s ; ETA = 0.63s Iter = 111 ; Dt = 0.0099s ; Time = 1.07s ; ETA = 0.60s Iter = 115 ; Dt = 0.0099s ; Time = 1.11s ; ETA = 0.59s Iter = 119 ; Dt = 0.0099s ; Time = 1.15s ; ETA = 0.55s Iter = 124 ; Dt = 0.0099s ; Time = 1.20s ; ETA = 0.53s Iter = 128 ; Dt = 0.0099s ; Time = 1.24s ; ETA = 0.51s Iter = 132 ; Dt = 0.0099s ; Time = 1.28s ; ETA = 0.47s Iter = 136 ; Dt = 0.0099s ; Time = 1.32s ; ETA = 0.45s Iter = 140 ; Dt = 0.0099s ; Time = 1.35s ; ETA = 0.42s Iter = 144 ; Dt = 0.0099s ; Time = 1.39s ; ETA = 0.41s Iter = 148 ; Dt = 0.0099s ; Time = 1.43s ; ETA = 0.38s Iter = 152 ; Dt = 0.0100s ; Time = 1.47s ; ETA = 0.35s Iter = 157 ; Dt = 0.0100s ; Time = 1.52s ; ETA = 0.31s Iter = 161 ; Dt = 0.0100s ; Time = 1.56s ; ETA = 0.29s Iter = 165 ; Dt = 0.0100s ; Time = 1.60s ; ETA = 0.26s Iter = 169 ; Dt = 0.0100s ; Time = 1.64s ; ETA = 0.24s Iter = 173 ; Dt = 0.0101s ; Time = 1.68s ; ETA = 0.21s Iter = 177 ; Dt = 0.0101s ; Time = 1.72s ; ETA = 0.18s Iter = 181 ; Dt = 0.0101s ; Time = 1.77s ; ETA = 0.16s Iter = 185 ; Dt = 0.0101s ; Time = 1.81s ; ETA = 0.13s Iter = 189 ; Dt = 0.0102s ; Time = 1.85s ; ETA = 0.10s Iter = 193 ; Dt = 0.0102s ; Time = 1.89s ; ETA = 0.08s Iter = 197 ; Dt = 0.0102s ; Time = 1.93s ; ETA = 0.05s Iter = 201 ; Dt = 0.0103s ; Time = 1.97s ; ETA = 0.02s Iter = 205 ; Dt = 0.0060s ; Time = 2.01s ; ETA = 0.00s =================================================================== | END | =================================================================== Problem.solve() completed in 1.3431456089019775s (wall time) SOLVING CASE 1 =================================================================== | INITIALIZATION | =================================================================== Problem size: Nodes=463, Layers=10, Triangles=848, Iter = 0 ; Dt = 0.0000s ; Time = 0.00s ; ETA = 0.00s =================================================================== | TIME LOOP | =================================================================== Iter = 5 ; Dt = 0.0093s ; Time = 0.05s ; ETA = 4.07s Iter = 9 ; Dt = 0.0093s ; Time = 0.08s ; ETA = 3.91s Iter = 14 ; Dt = 0.0093s ; Time = 0.13s ; ETA = 3.88s Iter = 18 ; Dt = 0.0093s ; Time = 0.17s ; ETA = 3.77s Iter = 22 ; Dt = 0.0093s ; Time = 0.21s ; ETA = 3.72s Iter = 27 ; Dt = 0.0093s ; Time = 0.25s ; ETA = 3.61s Iter = 31 ; Dt = 0.0093s ; Time = 0.29s ; ETA = 3.57s Iter = 36 ; Dt = 0.0093s ; Time = 0.34s ; ETA = 3.44s Iter = 40 ; Dt = 0.0093s ; Time = 0.37s ; ETA = 3.42s Iter = 44 ; Dt = 0.0093s ; Time = 0.41s ; ETA = 3.38s Iter = 49 ; Dt = 0.0094s ; Time = 0.46s ; ETA = 3.30s Iter = 53 ; Dt = 0.0094s ; Time = 0.49s ; ETA = 3.21s Iter = 58 ; Dt = 0.0087s ; Time = 0.54s ; ETA = 3.39s Iter = 62 ; Dt = 0.0082s ; Time = 0.57s ; ETA = 3.54s Iter = 68 ; Dt = 0.0077s ; Time = 0.62s ; ETA = 3.63s Iter = 73 ; Dt = 0.0074s ; Time = 0.66s ; ETA = 3.70s Iter = 79 ; Dt = 0.0071s ; Time = 0.70s ; ETA = 3.72s Iter = 85 ; Dt = 0.0069s ; Time = 0.74s ; ETA = 3.72s Iter = 91 ; Dt = 0.0068s ; Time = 0.78s ; ETA = 3.58s Iter = 97 ; Dt = 0.0068s ; Time = 0.83s ; ETA = 3.59s Iter = 102 ; Dt = 0.0070s ; Time = 0.86s ; ETA = 3.28s Iter = 108 ; Dt = 0.0081s ; Time = 0.91s ; ETA = 2.75s Iter = 113 ; Dt = 0.0078s ; Time = 0.95s ; ETA = 2.79s Iter = 118 ; Dt = 0.0080s ; Time = 0.99s ; ETA = 2.64s Iter = 123 ; Dt = 0.0072s ; Time = 1.02s ; ETA = 2.83s Iter = 130 ; Dt = 0.0070s ; Time = 1.07s ; ETA = 2.77s Iter = 136 ; Dt = 0.0070s ; Time = 1.11s ; ETA = 2.63s Iter = 143 ; Dt = 0.0054s ; Time = 1.15s ; ETA = 3.26s Iter = 149 ; Dt = 0.0069s ; Time = 1.19s ; ETA = 2.43s Iter = 155 ; Dt = 0.0074s ; Time = 1.23s ; ETA = 2.17s Iter = 161 ; Dt = 0.0078s ; Time = 1.28s ; ETA = 1.95s Iter = 166 ; Dt = 0.0077s ; Time = 1.32s ; ETA = 1.88s Iter = 171 ; Dt = 0.0077s ; Time = 1.35s ; ETA = 1.74s Iter = 177 ; Dt = 0.0069s ; Time = 1.40s ; ETA = 1.84s Iter = 183 ; Dt = 0.0069s ; Time = 1.44s ; ETA = 1.73s Iter = 189 ; Dt = 0.0068s ; Time = 1.48s ; ETA = 1.73s Iter = 195 ; Dt = 0.0076s ; Time = 1.52s ; ETA = 1.31s Iter = 200 ; Dt = 0.0077s ; Time = 1.56s ; ETA = 1.19s Iter = 206 ; Dt = 0.0080s ; Time = 1.60s ; ETA = 1.03s Iter = 211 ; Dt = 0.0082s ; Time = 1.64s ; ETA = 0.93s Iter = 216 ; Dt = 0.0079s ; Time = 1.68s ; ETA = 0.86s Iter = 221 ; Dt = 0.0078s ; Time = 1.72s ; ETA = 0.75s Iter = 227 ; Dt = 0.0071s ; Time = 1.76s ; ETA = 0.72s Iter = 233 ; Dt = 0.0070s ; Time = 1.80s ; ETA = 0.60s Iter = 239 ; Dt = 0.0069s ; Time = 1.85s ; ETA = 0.48s Iter = 245 ; Dt = 0.0083s ; Time = 1.89s ; ETA = 0.30s Iter = 250 ; Dt = 0.0083s ; Time = 1.93s ; ETA = 0.19s Iter = 255 ; Dt = 0.0082s ; Time = 1.97s ; ETA = 0.09s Iter = 260 ; Dt = 0.0028s ; Time = 2.01s ; ETA = 0.00s =================================================================== | END | =================================================================== Problem.solve() completed in 5.349927663803101s (wall time) SOLVING CASE 2 =================================================================== | INITIALIZATION | =================================================================== Problem size: Nodes=1054, Layers=12, Triangles=1990, Iter = 0 ; Dt = 0.0000s ; Time = 0.00s ; ETA = 0.00s =================================================================== | TIME LOOP | =================================================================== Iter = 7 ; Dt = 0.0061s ; Time = 0.04s ; ETA = 5.63s Iter = 14 ; Dt = 0.0061s ; Time = 0.09s ; ETA = 4.91s Iter = 21 ; Dt = 0.0061s ; Time = 0.13s ; ETA = 4.81s Iter = 27 ; Dt = 0.0061s ; Time = 0.17s ; ETA = 4.74s Iter = 34 ; Dt = 0.0061s ; Time = 0.21s ; ETA = 4.61s Iter = 41 ; Dt = 0.0062s ; Time = 0.25s ; ETA = 4.40s Iter = 47 ; Dt = 0.0062s ; Time = 0.29s ; ETA = 4.37s Iter = 54 ; Dt = 0.0062s ; Time = 0.33s ; ETA = 4.23s Iter = 61 ; Dt = 0.0062s ; Time = 0.37s ; ETA = 4.21s Iter = 67 ; Dt = 0.0062s ; Time = 0.41s ; ETA = 4.06s Iter = 74 ; Dt = 0.0062s ; Time = 0.46s ; ETA = 3.97s Iter = 80 ; Dt = 0.0062s ; Time = 0.49s ; ETA = 3.84s Iter = 87 ; Dt = 0.0062s ; Time = 0.54s ; ETA = 3.76s Iter = 93 ; Dt = 0.0063s ; Time = 0.57s ; ETA = 3.75s Iter = 100 ; Dt = 0.0063s ; Time = 0.62s ; ETA = 3.59s Iter = 107 ; Dt = 0.0063s ; Time = 0.66s ; ETA = 3.53s Iter = 113 ; Dt = 0.0063s ; Time = 0.70s ; ETA = 3.44s Iter = 119 ; Dt = 0.0063s ; Time = 0.74s ; ETA = 3.31s Iter = 126 ; Dt = 0.0063s ; Time = 0.78s ; ETA = 3.21s Iter = 132 ; Dt = 0.0063s ; Time = 0.82s ; ETA = 3.14s Iter = 139 ; Dt = 0.0064s ; Time = 0.86s ; ETA = 3.00s Iter = 145 ; Dt = 0.0064s ; Time = 0.90s ; ETA = 2.91s Iter = 152 ; Dt = 0.0064s ; Time = 0.95s ; ETA = 2.79s Iter = 158 ; Dt = 0.0064s ; Time = 0.98s ; ETA = 2.73s Iter = 165 ; Dt = 0.0064s ; Time = 1.03s ; ETA = 2.63s Iter = 171 ; Dt = 0.0064s ; Time = 1.07s ; ETA = 2.55s Iter = 177 ; Dt = 0.0064s ; Time = 1.11s ; ETA = 2.40s Iter = 184 ; Dt = 0.0064s ; Time = 1.15s ; ETA = 2.30s Iter = 190 ; Dt = 0.0064s ; Time = 1.19s ; ETA = 2.19s Iter = 197 ; Dt = 0.0064s ; Time = 1.23s ; ETA = 2.06s Iter = 203 ; Dt = 0.0064s ; Time = 1.27s ; ETA = 1.97s Iter = 210 ; Dt = 0.0064s ; Time = 1.32s ; ETA = 1.84s Iter = 216 ; Dt = 0.0064s ; Time = 1.35s ; ETA = 1.72s Iter = 222 ; Dt = 0.0064s ; Time = 1.39s ; ETA = 1.63s Iter = 229 ; Dt = 0.0064s ; Time = 1.44s ; ETA = 1.56s Iter = 235 ; Dt = 0.0064s ; Time = 1.48s ; ETA = 1.43s Iter = 242 ; Dt = 0.0064s ; Time = 1.52s ; ETA = 1.30s Iter = 248 ; Dt = 0.0064s ; Time = 1.56s ; ETA = 1.17s Iter = 254 ; Dt = 0.0064s ; Time = 1.60s ; ETA = 1.09s Iter = 261 ; Dt = 0.0064s ; Time = 1.64s ; ETA = 0.98s Iter = 267 ; Dt = 0.0064s ; Time = 1.68s ; ETA = 0.86s Iter = 273 ; Dt = 0.0065s ; Time = 1.72s ; ETA = 0.76s Iter = 280 ; Dt = 0.0065s ; Time = 1.76s ; ETA = 0.63s Iter = 286 ; Dt = 0.0065s ; Time = 1.80s ; ETA = 0.53s Iter = 292 ; Dt = 0.0065s ; Time = 1.84s ; ETA = 0.43s Iter = 299 ; Dt = 0.0065s ; Time = 1.89s ; ETA = 0.31s Iter = 305 ; Dt = 0.0065s ; Time = 1.93s ; ETA = 0.21s Iter = 311 ; Dt = 0.0065s ; Time = 1.97s ; ETA = 0.11s Iter = 318 ; Dt = 0.0006s ; Time = 2.01s ; ETA = 0.00s =================================================================== | END | =================================================================== Problem.solve() completed in 5.425991058349609s (wall time) SOLVING CASE 3 =================================================================== | INITIALIZATION | =================================================================== Problem size: Nodes=1054, Layers=12, Triangles=1990, Iter = 0 ; Dt = 0.0000s ; Time = 0.00s ; ETA = 0.00s =================================================================== | TIME LOOP | =================================================================== Iter = 7 ; Dt = 0.0061s ; Time = 0.04s ; ETA = 16.62s Iter = 14 ; Dt = 0.0061s ; Time = 0.09s ; ETA = 15.97s Iter = 21 ; Dt = 0.0061s ; Time = 0.13s ; ETA = 15.73s Iter = 27 ; Dt = 0.0061s ; Time = 0.16s ; ETA = 15.69s Iter = 34 ; Dt = 0.0061s ; Time = 0.21s ; ETA = 15.04s Iter = 41 ; Dt = 0.0061s ; Time = 0.25s ; ETA = 15.00s Iter = 47 ; Dt = 0.0061s ; Time = 0.29s ; ETA = 14.68s Iter = 54 ; Dt = 0.0061s ; Time = 0.33s ; ETA = 15.74s Iter = 61 ; Dt = 0.0061s ; Time = 0.37s ; ETA = 14.21s Iter = 68 ; Dt = 0.0061s ; Time = 0.41s ; ETA = 14.04s Iter = 74 ; Dt = 0.0061s ; Time = 0.45s ; ETA = 13.82s Iter = 81 ; Dt = 0.0061s ; Time = 0.49s ; ETA = 13.68s Iter = 88 ; Dt = 0.0056s ; Time = 0.53s ; ETA = 13.97s Iter = 96 ; Dt = 0.0053s ; Time = 0.58s ; ETA = 13.99s Iter = 104 ; Dt = 0.0049s ; Time = 0.62s ; ETA = 14.76s Iter = 112 ; Dt = 0.0049s ; Time = 0.66s ; ETA = 14.52s Iter = 121 ; Dt = 0.0049s ; Time = 0.70s ; ETA = 13.85s Iter = 129 ; Dt = 0.0049s ; Time = 0.74s ; ETA = 13.54s Iter = 137 ; Dt = 0.0051s ; Time = 0.78s ; ETA = 12.45s Iter = 145 ; Dt = 0.0051s ; Time = 0.82s ; ETA = 12.36s Iter = 153 ; Dt = 0.0053s ; Time = 0.86s ; ETA = 11.56s Iter = 161 ; Dt = 0.0052s ; Time = 0.90s ; ETA = 11.36s Iter = 169 ; Dt = 0.0048s ; Time = 0.94s ; ETA = 12.17s Iter = 178 ; Dt = 0.0048s ; Time = 0.99s ; ETA = 11.43s Iter = 187 ; Dt = 0.0047s ; Time = 1.03s ; ETA = 11.41s Iter = 195 ; Dt = 0.0049s ; Time = 1.07s ; ETA = 10.53s Iter = 204 ; Dt = 0.0047s ; Time = 1.11s ; ETA = 10.43s Iter = 213 ; Dt = 0.0050s ; Time = 1.15s ; ETA = 9.48s Iter = 221 ; Dt = 0.0052s ; Time = 1.19s ; ETA = 8.65s Iter = 229 ; Dt = 0.0055s ; Time = 1.23s ; ETA = 7.72s Iter = 236 ; Dt = 0.0058s ; Time = 1.27s ; ETA = 7.09s Iter = 244 ; Dt = 0.0042s ; Time = 1.31s ; ETA = 9.28s Iter = 252 ; Dt = 0.0053s ; Time = 1.36s ; ETA = 6.92s Iter = 260 ; Dt = 0.0035s ; Time = 1.39s ; ETA = 9.76s Iter = 269 ; Dt = 0.0043s ; Time = 1.44s ; ETA = 7.40s Iter = 277 ; Dt = 0.0052s ; Time = 1.47s ; ETA = 5.58s Iter = 285 ; Dt = 0.0055s ; Time = 1.52s ; ETA = 4.97s Iter = 293 ; Dt = 0.0055s ; Time = 1.56s ; ETA = 4.48s Iter = 301 ; Dt = 0.0052s ; Time = 1.60s ; ETA = 4.25s Iter = 309 ; Dt = 0.0053s ; Time = 1.64s ; ETA = 3.73s Iter = 318 ; Dt = 0.0049s ; Time = 1.68s ; ETA = 3.63s Iter = 326 ; Dt = 0.0046s ; Time = 1.72s ; ETA = 3.40s Iter = 335 ; Dt = 0.0055s ; Time = 1.77s ; ETA = 2.45s Iter = 343 ; Dt = 0.0055s ; Time = 1.80s ; ETA = 2.09s Iter = 351 ; Dt = 0.0048s ; Time = 1.84s ; ETA = 1.87s Iter = 360 ; Dt = 0.0048s ; Time = 1.89s ; ETA = 1.39s Iter = 369 ; Dt = 0.0051s ; Time = 1.93s ; ETA = 0.86s Iter = 377 ; Dt = 0.0052s ; Time = 1.97s ; ETA = 0.39s Iter = 385 ; Dt = 0.0020s ; Time = 2.01s ; ETA = 0.00s =================================================================== | END | =================================================================== Problem.solve() completed in 20.98289942741394s (wall time) SOLVING CASE 4 =================================================================== | INITIALIZATION | =================================================================== Problem size: Nodes=1848, Layers=14, Triangles=3538, Iter = 0 ; Dt = 0.0000s ; Time = 0.00s ; ETA = 0.00s =================================================================== | TIME LOOP | =================================================================== Iter = 10 ; Dt = 0.0045s ; Time = 0.05s ; ETA = 13.95s Iter = 19 ; Dt = 0.0045s ; Time = 0.09s ; ETA = 12.91s Iter = 28 ; Dt = 0.0045s ; Time = 0.13s ; ETA = 12.41s Iter = 37 ; Dt = 0.0046s ; Time = 0.17s ; ETA = 12.36s Iter = 46 ; Dt = 0.0046s ; Time = 0.21s ; ETA = 12.02s Iter = 55 ; Dt = 0.0046s ; Time = 0.25s ; ETA = 11.73s Iter = 63 ; Dt = 0.0046s ; Time = 0.29s ; ETA = 11.62s Iter = 72 ; Dt = 0.0046s ; Time = 0.33s ; ETA = 11.10s Iter = 81 ; Dt = 0.0046s ; Time = 0.37s ; ETA = 10.79s Iter = 90 ; Dt = 0.0046s ; Time = 0.41s ; ETA = 10.95s Iter = 99 ; Dt = 0.0046s ; Time = 0.45s ; ETA = 10.60s Iter = 108 ; Dt = 0.0046s ; Time = 0.49s ; ETA = 10.26s Iter = 117 ; Dt = 0.0046s ; Time = 0.53s ; ETA = 10.08s Iter = 126 ; Dt = 0.0046s ; Time = 0.58s ; ETA = 9.63s Iter = 135 ; Dt = 0.0046s ; Time = 0.62s ; ETA = 9.56s Iter = 143 ; Dt = 0.0046s ; Time = 0.66s ; ETA = 9.36s Iter = 152 ; Dt = 0.0047s ; Time = 0.70s ; ETA = 9.16s Iter = 161 ; Dt = 0.0047s ; Time = 0.74s ; ETA = 8.90s Iter = 170 ; Dt = 0.0047s ; Time = 0.78s ; ETA = 8.48s Iter = 179 ; Dt = 0.0047s ; Time = 0.82s ; ETA = 8.16s Iter = 187 ; Dt = 0.0047s ; Time = 0.86s ; ETA = 7.88s Iter = 196 ; Dt = 0.0047s ; Time = 0.90s ; ETA = 7.71s Iter = 205 ; Dt = 0.0047s ; Time = 0.95s ; ETA = 7.47s Iter = 213 ; Dt = 0.0047s ; Time = 0.98s ; ETA = 7.21s Iter = 222 ; Dt = 0.0047s ; Time = 1.02s ; ETA = 6.97s Iter = 231 ; Dt = 0.0047s ; Time = 1.07s ; ETA = 6.57s Iter = 240 ; Dt = 0.0047s ; Time = 1.11s ; ETA = 6.43s Iter = 248 ; Dt = 0.0047s ; Time = 1.15s ; ETA = 6.02s Iter = 257 ; Dt = 0.0047s ; Time = 1.19s ; ETA = 5.68s Iter = 266 ; Dt = 0.0047s ; Time = 1.23s ; ETA = 5.44s Iter = 275 ; Dt = 0.0047s ; Time = 1.27s ; ETA = 5.11s Iter = 283 ; Dt = 0.0047s ; Time = 1.31s ; ETA = 4.97s Iter = 292 ; Dt = 0.0047s ; Time = 1.35s ; ETA = 4.61s Iter = 301 ; Dt = 0.0047s ; Time = 1.40s ; ETA = 4.32s Iter = 310 ; Dt = 0.0047s ; Time = 1.44s ; ETA = 4.07s Iter = 318 ; Dt = 0.0047s ; Time = 1.47s ; ETA = 3.78s Iter = 327 ; Dt = 0.0047s ; Time = 1.52s ; ETA = 3.47s Iter = 336 ; Dt = 0.0047s ; Time = 1.56s ; ETA = 3.21s Iter = 344 ; Dt = 0.0047s ; Time = 1.60s ; ETA = 2.91s Iter = 353 ; Dt = 0.0047s ; Time = 1.64s ; ETA = 2.61s Iter = 362 ; Dt = 0.0047s ; Time = 1.68s ; ETA = 2.34s Iter = 370 ; Dt = 0.0047s ; Time = 1.72s ; ETA = 2.02s Iter = 379 ; Dt = 0.0047s ; Time = 1.76s ; ETA = 1.76s Iter = 388 ; Dt = 0.0047s ; Time = 1.81s ; ETA = 1.45s Iter = 396 ; Dt = 0.0048s ; Time = 1.84s ; ETA = 1.16s Iter = 405 ; Dt = 0.0048s ; Time = 1.89s ; ETA = 0.86s Iter = 414 ; Dt = 0.0048s ; Time = 1.93s ; ETA = 0.54s Iter = 422 ; Dt = 0.0048s ; Time = 1.97s ; ETA = 0.28s Iter = 431 ; Dt = 0.0007s ; Time = 2.01s ; ETA = 0.00s =================================================================== | END | =================================================================== Problem.solve() completed in 14.246922492980957s (wall time) SOLVING CASE 5 =================================================================== | INITIALIZATION | =================================================================== Problem size: Nodes=1848, Layers=14, Triangles=3538, Iter = 0 ; Dt = 0.0000s ; Time = 0.00s ; ETA = 0.00s =================================================================== | TIME LOOP | =================================================================== Iter = 10 ; Dt = 0.0045s ; Time = 0.05s ; ETA = 46.18s Iter = 19 ; Dt = 0.0045s ; Time = 0.09s ; ETA = 45.14s Iter = 28 ; Dt = 0.0045s ; Time = 0.13s ; ETA = 42.68s Iter = 37 ; Dt = 0.0045s ; Time = 0.17s ; ETA = 41.48s Iter = 46 ; Dt = 0.0045s ; Time = 0.21s ; ETA = 42.27s Iter = 55 ; Dt = 0.0045s ; Time = 0.25s ; ETA = 41.20s Iter = 64 ; Dt = 0.0045s ; Time = 0.29s ; ETA = 40.36s Iter = 73 ; Dt = 0.0045s ; Time = 0.33s ; ETA = 39.57s Iter = 82 ; Dt = 0.0045s ; Time = 0.37s ; ETA = 42.07s Iter = 91 ; Dt = 0.0045s ; Time = 0.41s ; ETA = 39.03s Iter = 100 ; Dt = 0.0041s ; Time = 0.45s ; ETA = 42.32s Iter = 111 ; Dt = 0.0038s ; Time = 0.49s ; ETA = 43.02s Iter = 121 ; Dt = 0.0037s ; Time = 0.53s ; ETA = 43.11s Iter = 132 ; Dt = 0.0037s ; Time = 0.57s ; ETA = 42.80s Iter = 143 ; Dt = 0.0037s ; Time = 0.62s ; ETA = 42.86s Iter = 154 ; Dt = 0.0040s ; Time = 0.66s ; ETA = 37.90s Iter = 164 ; Dt = 0.0042s ; Time = 0.70s ; ETA = 35.37s Iter = 173 ; Dt = 0.0041s ; Time = 0.74s ; ETA = 35.96s Iter = 183 ; Dt = 0.0041s ; Time = 0.78s ; ETA = 33.27s Iter = 193 ; Dt = 0.0039s ; Time = 0.82s ; ETA = 34.60s Iter = 204 ; Dt = 0.0037s ; Time = 0.86s ; ETA = 34.85s Iter = 214 ; Dt = 0.0038s ; Time = 0.90s ; ETA = 33.78s Iter = 226 ; Dt = 0.0037s ; Time = 0.95s ; ETA = 32.38s Iter = 236 ; Dt = 0.0039s ; Time = 0.98s ; ETA = 30.51s Iter = 247 ; Dt = 0.0038s ; Time = 1.03s ; ETA = 29.20s Iter = 257 ; Dt = 0.0041s ; Time = 1.06s ; ETA = 26.70s Iter = 267 ; Dt = 0.0042s ; Time = 1.11s ; ETA = 25.39s Iter = 277 ; Dt = 0.0043s ; Time = 1.15s ; ETA = 23.55s Iter = 287 ; Dt = 0.0040s ; Time = 1.19s ; ETA = 23.51s Iter = 298 ; Dt = 0.0037s ; Time = 1.23s ; ETA = 24.06s Iter = 310 ; Dt = 0.0038s ; Time = 1.27s ; ETA = 21.92s Iter = 321 ; Dt = 0.0037s ; Time = 1.31s ; ETA = 21.95s Iter = 333 ; Dt = 0.0026s ; Time = 1.35s ; ETA = 28.80s Iter = 345 ; Dt = 0.0043s ; Time = 1.39s ; ETA = 16.66s Iter = 355 ; Dt = 0.0039s ; Time = 1.43s ; ETA = 17.32s Iter = 368 ; Dt = 0.0035s ; Time = 1.48s ; ETA = 17.52s Iter = 378 ; Dt = 0.0041s ; Time = 1.52s ; ETA = 13.89s Iter = 390 ; Dt = 0.0034s ; Time = 1.56s ; ETA = 15.47s Iter = 401 ; Dt = 0.0040s ; Time = 1.60s ; ETA = 11.67s Iter = 413 ; Dt = 0.0031s ; Time = 1.64s ; ETA = 13.93s Iter = 424 ; Dt = 0.0040s ; Time = 1.68s ; ETA = 9.46s Iter = 434 ; Dt = 0.0043s ; Time = 1.72s ; ETA = 7.50s Iter = 447 ; Dt = 0.0043s ; Time = 1.76s ; ETA = 6.67s Iter = 457 ; Dt = 0.0033s ; Time = 1.80s ; ETA = 7.13s Iter = 470 ; Dt = 0.0042s ; Time = 1.84s ; ETA = 4.52s Iter = 482 ; Dt = 0.0038s ; Time = 1.88s ; ETA = 3.79s Iter = 494 ; Dt = 0.0021s ; Time = 1.93s ; ETA = 4.51s Iter = 506 ; Dt = 0.0042s ; Time = 1.97s ; ETA = 1.14s Iter = 516 ; Dt = 0.0036s ; Time = 2.01s ; ETA = 0.00s =================================================================== | END | =================================================================== Problem.solve() completed in 58.85429286956787s (wall time) SOLVING CASE 6 =================================================================== | INITIALIZATION | =================================================================== Problem size: Nodes=2914, Layers=16, Triangles=5630, Iter = 0 ; Dt = 0.0000s ; Time = 0.00s ; ETA = 0.00s =================================================================== | TIME LOOP | =================================================================== Iter = 12 ; Dt = 0.0036s ; Time = 0.04s ; ETA = 32.00s Iter = 23 ; Dt = 0.0036s ; Time = 0.08s ; ETA = 30.30s Iter = 34 ; Dt = 0.0036s ; Time = 0.12s ; ETA = 29.57s Iter = 46 ; Dt = 0.0036s ; Time = 0.17s ; ETA = 28.68s Iter = 57 ; Dt = 0.0036s ; Time = 0.21s ; ETA = 28.10s Iter = 68 ; Dt = 0.0036s ; Time = 0.25s ; ETA = 27.63s Iter = 80 ; Dt = 0.0036s ; Time = 0.29s ; ETA = 27.56s Iter = 91 ; Dt = 0.0036s ; Time = 0.33s ; ETA = 26.26s Iter = 102 ; Dt = 0.0036s ; Time = 0.37s ; ETA = 25.78s Iter = 113 ; Dt = 0.0036s ; Time = 0.41s ; ETA = 25.50s Iter = 125 ; Dt = 0.0037s ; Time = 0.45s ; ETA = 24.75s Iter = 136 ; Dt = 0.0037s ; Time = 0.49s ; ETA = 24.46s Iter = 147 ; Dt = 0.0037s ; Time = 0.53s ; ETA = 23.71s Iter = 158 ; Dt = 0.0037s ; Time = 0.57s ; ETA = 23.27s Iter = 169 ; Dt = 0.0037s ; Time = 0.61s ; ETA = 23.14s Iter = 180 ; Dt = 0.0037s ; Time = 0.66s ; ETA = 21.76s Iter = 192 ; Dt = 0.0037s ; Time = 0.70s ; ETA = 21.40s Iter = 203 ; Dt = 0.0037s ; Time = 0.74s ; ETA = 21.19s Iter = 214 ; Dt = 0.0037s ; Time = 0.78s ; ETA = 20.17s Iter = 225 ; Dt = 0.0037s ; Time = 0.82s ; ETA = 19.24s Iter = 236 ; Dt = 0.0037s ; Time = 0.86s ; ETA = 18.31s Iter = 247 ; Dt = 0.0037s ; Time = 0.90s ; ETA = 18.64s Iter = 258 ; Dt = 0.0037s ; Time = 0.94s ; ETA = 18.22s Iter = 269 ; Dt = 0.0037s ; Time = 0.98s ; ETA = 17.49s Iter = 280 ; Dt = 0.0037s ; Time = 1.03s ; ETA = 15.61s Iter = 291 ; Dt = 0.0037s ; Time = 1.07s ; ETA = 15.07s Iter = 302 ; Dt = 0.0037s ; Time = 1.11s ; ETA = 15.35s Iter = 313 ; Dt = 0.0037s ; Time = 1.15s ; ETA = 15.03s Iter = 324 ; Dt = 0.0037s ; Time = 1.19s ; ETA = 14.23s Iter = 335 ; Dt = 0.0037s ; Time = 1.23s ; ETA = 13.51s Iter = 346 ; Dt = 0.0037s ; Time = 1.27s ; ETA = 11.80s Iter = 357 ; Dt = 0.0037s ; Time = 1.31s ; ETA = 11.27s Iter = 368 ; Dt = 0.0037s ; Time = 1.35s ; ETA = 11.22s Iter = 379 ; Dt = 0.0037s ; Time = 1.39s ; ETA = 10.90s Iter = 391 ; Dt = 0.0037s ; Time = 1.44s ; ETA = 10.01s Iter = 402 ; Dt = 0.0037s ; Time = 1.48s ; ETA = 9.49s Iter = 413 ; Dt = 0.0037s ; Time = 1.52s ; ETA = 8.63s Iter = 424 ; Dt = 0.0037s ; Time = 1.56s ; ETA = 7.94s Iter = 435 ; Dt = 0.0037s ; Time = 1.60s ; ETA = 6.61s Iter = 446 ; Dt = 0.0037s ; Time = 1.64s ; ETA = 5.88s Iter = 457 ; Dt = 0.0037s ; Time = 1.68s ; ETA = 5.29s Iter = 468 ; Dt = 0.0037s ; Time = 1.72s ; ETA = 4.82s Iter = 479 ; Dt = 0.0037s ; Time = 1.76s ; ETA = 4.18s Iter = 489 ; Dt = 0.0037s ; Time = 1.80s ; ETA = 3.54s Iter = 500 ; Dt = 0.0037s ; Time = 1.84s ; ETA = 2.74s Iter = 511 ; Dt = 0.0038s ; Time = 1.88s ; ETA = 2.04s Iter = 522 ; Dt = 0.0038s ; Time = 1.93s ; ETA = 1.39s Iter = 533 ; Dt = 0.0038s ; Time = 1.97s ; ETA = 0.63s Iter = 544 ; Dt = 0.0017s ; Time = 2.01s ; ETA = 0.00s =================================================================== | END | =================================================================== Problem.solve() completed in 33.53608989715576s (wall time) SOLVING CASE 7 =================================================================== | INITIALIZATION | =================================================================== Problem size: Nodes=2914, Layers=16, Triangles=5630, Iter = 0 ; Dt = 0.0000s ; Time = 0.00s ; ETA = 0.00s =================================================================== | TIME LOOP | =================================================================== Iter = 12 ; Dt = 0.0036s ; Time = 0.04s ; ETA = 118.39s Iter = 23 ; Dt = 0.0036s ; Time = 0.08s ; ETA = 110.71s Iter = 35 ; Dt = 0.0036s ; Time = 0.13s ; ETA = 111.10s Iter = 46 ; Dt = 0.0036s ; Time = 0.17s ; ETA = 108.06s Iter = 57 ; Dt = 0.0036s ; Time = 0.21s ; ETA = 111.38s Iter = 69 ; Dt = 0.0036s ; Time = 0.25s ; ETA = 105.40s Iter = 80 ; Dt = 0.0036s ; Time = 0.29s ; ETA = 101.64s Iter = 91 ; Dt = 0.0036s ; Time = 0.33s ; ETA = 98.05s Iter = 103 ; Dt = 0.0036s ; Time = 0.37s ; ETA = 92.92s Iter = 114 ; Dt = 0.0033s ; Time = 0.41s ; ETA = 105.07s Iter = 127 ; Dt = 0.0030s ; Time = 0.45s ; ETA = 107.86s Iter = 141 ; Dt = 0.0032s ; Time = 0.49s ; ETA = 98.11s Iter = 154 ; Dt = 0.0032s ; Time = 0.54s ; ETA = 104.81s Iter = 166 ; Dt = 0.0033s ; Time = 0.57s ; ETA = 95.74s Iter = 178 ; Dt = 0.0034s ; Time = 0.62s ; ETA = 87.35s Iter = 190 ; Dt = 0.0036s ; Time = 0.66s ; ETA = 90.90s Iter = 201 ; Dt = 0.0033s ; Time = 0.70s ; ETA = 89.83s Iter = 214 ; Dt = 0.0032s ; Time = 0.74s ; ETA = 89.69s Iter = 226 ; Dt = 0.0032s ; Time = 0.78s ; ETA = 88.24s Iter = 239 ; Dt = 0.0031s ; Time = 0.82s ; ETA = 90.87s Iter = 252 ; Dt = 0.0030s ; Time = 0.86s ; ETA = 87.43s Iter = 265 ; Dt = 0.0033s ; Time = 0.90s ; ETA = 73.91s Iter = 277 ; Dt = 0.0034s ; Time = 0.94s ; ETA = 68.09s Iter = 289 ; Dt = 0.0034s ; Time = 0.98s ; ETA = 69.81s Iter = 301 ; Dt = 0.0035s ; Time = 1.02s ; ETA = 64.43s Iter = 313 ; Dt = 0.0034s ; Time = 1.07s ; ETA = 62.92s Iter = 326 ; Dt = 0.0034s ; Time = 1.11s ; ETA = 59.03s Iter = 339 ; Dt = 0.0031s ; Time = 1.15s ; ETA = 61.38s Iter = 353 ; Dt = 0.0030s ; Time = 1.19s ; ETA = 59.35s Iter = 367 ; Dt = 0.0029s ; Time = 1.23s ; ETA = 60.31s Iter = 381 ; Dt = 0.0027s ; Time = 1.27s ; ETA = 62.98s Iter = 396 ; Dt = 0.0017s ; Time = 1.31s ; ETA = 91.19s Iter = 409 ; Dt = 0.0035s ; Time = 1.35s ; ETA = 42.73s Iter = 424 ; Dt = 0.0032s ; Time = 1.39s ; ETA = 46.22s Iter = 438 ; Dt = 0.0034s ; Time = 1.43s ; ETA = 40.48s Iter = 452 ; Dt = 0.0025s ; Time = 1.47s ; ETA = 51.59s Iter = 466 ; Dt = 0.0030s ; Time = 1.52s ; ETA = 38.65s Iter = 481 ; Dt = 0.0023s ; Time = 1.56s ; ETA = 45.16s Iter = 496 ; Dt = 0.0030s ; Time = 1.60s ; ETA = 34.09s Iter = 511 ; Dt = 0.0025s ; Time = 1.64s ; ETA = 33.16s Iter = 527 ; Dt = 0.0036s ; Time = 1.68s ; ETA = 21.41s Iter = 543 ; Dt = 0.0021s ; Time = 1.72s ; ETA = 30.59s Iter = 558 ; Dt = 0.0016s ; Time = 1.76s ; ETA = 33.73s Iter = 574 ; Dt = 0.0024s ; Time = 1.80s ; ETA = 20.50s Iter = 589 ; Dt = 0.0024s ; Time = 1.84s ; ETA = 15.10s Iter = 608 ; Dt = 0.0020s ; Time = 1.88s ; ETA = 14.27s Iter = 625 ; Dt = 0.0020s ; Time = 1.93s ; ETA = 9.60s Iter = 642 ; Dt = 0.0018s ; Time = 1.97s ; ETA = 5.43s Iter = 656 ; Dt = 0.0017s ; Time = 2.01s ; ETA = 0.00s =================================================================== | END | =================================================================== Problem.solve() completed in 149.6731457710266s (wall time) SOLVING CASE 8 =================================================================== | INITIALIZATION | =================================================================== Problem size: Nodes=4155, Layers=18, Triangles=8072, Iter = 0 ; Dt = 0.0000s ; Time = 0.00s ; ETA = 0.00s =================================================================== | TIME LOOP | =================================================================== Iter = 14 ; Dt = 0.0030s ; Time = 0.04s ; ETA = 58.38s Iter = 28 ; Dt = 0.0030s ; Time = 0.08s ; ETA = 56.06s Iter = 41 ; Dt = 0.0030s ; Time = 0.12s ; ETA = 54.77s Iter = 55 ; Dt = 0.0030s ; Time = 0.17s ; ETA = 54.16s Iter = 69 ; Dt = 0.0030s ; Time = 0.21s ; ETA = 53.07s Iter = 82 ; Dt = 0.0030s ; Time = 0.25s ; ETA = 51.15s Iter = 96 ; Dt = 0.0030s ; Time = 0.29s ; ETA = 51.23s Iter = 109 ; Dt = 0.0030s ; Time = 0.33s ; ETA = 49.87s Iter = 123 ; Dt = 0.0030s ; Time = 0.37s ; ETA = 48.40s Iter = 137 ; Dt = 0.0030s ; Time = 0.41s ; ETA = 48.01s Iter = 150 ; Dt = 0.0030s ; Time = 0.45s ; ETA = 47.11s Iter = 164 ; Dt = 0.0030s ; Time = 0.49s ; ETA = 44.92s Iter = 177 ; Dt = 0.0030s ; Time = 0.53s ; ETA = 44.39s Iter = 191 ; Dt = 0.0030s ; Time = 0.58s ; ETA = 43.66s Iter = 204 ; Dt = 0.0030s ; Time = 0.62s ; ETA = 41.75s Iter = 217 ; Dt = 0.0031s ; Time = 0.66s ; ETA = 41.43s Iter = 231 ; Dt = 0.0031s ; Time = 0.70s ; ETA = 39.25s Iter = 244 ; Dt = 0.0031s ; Time = 0.74s ; ETA = 38.59s Iter = 258 ; Dt = 0.0031s ; Time = 0.78s ; ETA = 37.33s Iter = 271 ; Dt = 0.0031s ; Time = 0.82s ; ETA = 36.32s Iter = 284 ; Dt = 0.0031s ; Time = 0.86s ; ETA = 34.72s Iter = 298 ; Dt = 0.0031s ; Time = 0.90s ; ETA = 33.58s Iter = 311 ; Dt = 0.0031s ; Time = 0.94s ; ETA = 32.95s Iter = 324 ; Dt = 0.0031s ; Time = 0.98s ; ETA = 31.04s Iter = 338 ; Dt = 0.0031s ; Time = 1.03s ; ETA = 29.84s Iter = 351 ; Dt = 0.0031s ; Time = 1.07s ; ETA = 28.75s Iter = 364 ; Dt = 0.0031s ; Time = 1.11s ; ETA = 28.44s Iter = 378 ; Dt = 0.0031s ; Time = 1.15s ; ETA = 27.06s Iter = 391 ; Dt = 0.0031s ; Time = 1.19s ; ETA = 25.45s Iter = 405 ; Dt = 0.0031s ; Time = 1.23s ; ETA = 24.44s Iter = 418 ; Dt = 0.0031s ; Time = 1.27s ; ETA = 23.04s Iter = 431 ; Dt = 0.0031s ; Time = 1.31s ; ETA = 21.63s Iter = 445 ; Dt = 0.0031s ; Time = 1.35s ; ETA = 20.62s Iter = 458 ; Dt = 0.0031s ; Time = 1.39s ; ETA = 19.28s Iter = 472 ; Dt = 0.0031s ; Time = 1.44s ; ETA = 17.97s Iter = 485 ; Dt = 0.0031s ; Time = 1.48s ; ETA = 16.80s Iter = 498 ; Dt = 0.0031s ; Time = 1.52s ; ETA = 15.59s Iter = 512 ; Dt = 0.0031s ; Time = 1.56s ; ETA = 14.18s Iter = 525 ; Dt = 0.0031s ; Time = 1.60s ; ETA = 12.92s Iter = 538 ; Dt = 0.0031s ; Time = 1.64s ; ETA = 11.48s Iter = 552 ; Dt = 0.0031s ; Time = 1.68s ; ETA = 10.29s Iter = 565 ; Dt = 0.0031s ; Time = 1.72s ; ETA = 8.93s Iter = 578 ; Dt = 0.0031s ; Time = 1.76s ; ETA = 7.65s Iter = 591 ; Dt = 0.0031s ; Time = 1.80s ; ETA = 6.38s Iter = 605 ; Dt = 0.0031s ; Time = 1.85s ; ETA = 5.23s Iter = 618 ; Dt = 0.0031s ; Time = 1.89s ; ETA = 3.80s Iter = 631 ; Dt = 0.0031s ; Time = 1.93s ; ETA = 2.50s Iter = 644 ; Dt = 0.0031s ; Time = 1.97s ; ETA = 1.27s Iter = 657 ; Dt = 0.0028s ; Time = 2.01s ; ETA = 0.00s =================================================================== | END | =================================================================== Problem.solve() completed in 62.273173332214355s (wall time) SOLVING CASE 9 =================================================================== | INITIALIZATION | =================================================================== Problem size: Nodes=4155, Layers=18, Triangles=8072, Iter = 0 ; Dt = 0.0000s ; Time = 0.00s ; ETA = 0.00s =================================================================== | TIME LOOP | =================================================================== Iter = 14 ; Dt = 0.0030s ; Time = 0.04s ; ETA = 221.60s Iter = 28 ; Dt = 0.0030s ; Time = 0.08s ; ETA = 217.25s Iter = 41 ; Dt = 0.0030s ; Time = 0.12s ; ETA = 215.09s Iter = 55 ; Dt = 0.0030s ; Time = 0.16s ; ETA = 212.14s Iter = 69 ; Dt = 0.0030s ; Time = 0.21s ; ETA = 208.36s Iter = 82 ; Dt = 0.0030s ; Time = 0.25s ; ETA = 204.83s Iter = 96 ; Dt = 0.0030s ; Time = 0.29s ; ETA = 207.64s Iter = 110 ; Dt = 0.0030s ; Time = 0.33s ; ETA = 205.21s Iter = 123 ; Dt = 0.0030s ; Time = 0.37s ; ETA = 202.61s Iter = 138 ; Dt = 0.0028s ; Time = 0.41s ; ETA = 201.13s Iter = 153 ; Dt = 0.0028s ; Time = 0.45s ; ETA = 218.65s Iter = 168 ; Dt = 0.0028s ; Time = 0.49s ; ETA = 191.92s Iter = 182 ; Dt = 0.0030s ; Time = 0.53s ; ETA = 180.46s Iter = 196 ; Dt = 0.0030s ; Time = 0.58s ; ETA = 172.72s Iter = 209 ; Dt = 0.0030s ; Time = 0.61s ; ETA = 170.25s Iter = 223 ; Dt = 0.0030s ; Time = 0.66s ; ETA = 167.58s Iter = 237 ; Dt = 0.0027s ; Time = 0.70s ; ETA = 180.39s Iter = 253 ; Dt = 0.0026s ; Time = 0.74s ; ETA = 185.02s Iter = 267 ; Dt = 0.0029s ; Time = 0.78s ; ETA = 159.11s Iter = 282 ; Dt = 0.0028s ; Time = 0.82s ; ETA = 156.77s Iter = 297 ; Dt = 0.0029s ; Time = 0.86s ; ETA = 158.78s Iter = 311 ; Dt = 0.0030s ; Time = 0.90s ; ETA = 149.79s Iter = 326 ; Dt = 0.0024s ; Time = 0.94s ; ETA = 169.98s Iter = 340 ; Dt = 0.0030s ; Time = 0.99s ; ETA = 129.38s Iter = 354 ; Dt = 0.0028s ; Time = 1.02s ; ETA = 132.97s Iter = 370 ; Dt = 0.0029s ; Time = 1.07s ; ETA = 123.72s Iter = 385 ; Dt = 0.0021s ; Time = 1.11s ; ETA = 166.46s Iter = 401 ; Dt = 0.0029s ; Time = 1.15s ; ETA = 115.10s Iter = 419 ; Dt = 0.0021s ; Time = 1.19s ; ETA = 153.36s Iter = 434 ; Dt = 0.0021s ; Time = 1.23s ; ETA = 146.23s Iter = 451 ; Dt = 0.0030s ; Time = 1.27s ; ETA = 103.05s Iter = 468 ; Dt = 0.0021s ; Time = 1.31s ; ETA = 124.63s Iter = 486 ; Dt = 0.0023s ; Time = 1.35s ; ETA = 111.50s Iter = 504 ; Dt = 0.0022s ; Time = 1.39s ; ETA = 119.46s Iter = 523 ; Dt = 0.0016s ; Time = 1.43s ; ETA = 141.62s Iter = 541 ; Dt = 0.0021s ; Time = 1.47s ; ETA = 96.67s Iter = 561 ; Dt = 0.0020s ; Time = 1.52s ; ETA = 93.68s Iter = 580 ; Dt = 0.0030s ; Time = 1.56s ; ETA = 57.43s Iter = 599 ; Dt = 0.0030s ; Time = 1.60s ; ETA = 52.29s Iter = 619 ; Dt = 0.0016s ; Time = 1.64s ; ETA = 86.33s Iter = 636 ; Dt = 0.0029s ; Time = 1.68s ; ETA = 43.04s Iter = 658 ; Dt = 0.0020s ; Time = 1.72s ; ETA = 55.41s Iter = 680 ; Dt = 0.0018s ; Time = 1.76s ; ETA = 60.19s Iter = 700 ; Dt = 0.0013s ; Time = 1.80s ; ETA = 60.48s Iter = 723 ; Dt = 0.0029s ; Time = 1.85s ; ETA = 21.51s Iter = 744 ; Dt = 0.0014s ; Time = 1.88s ; ETA = 32.83s Iter = 764 ; Dt = 0.0017s ; Time = 1.92s ; ETA = 18.38s Iter = 788 ; Dt = 0.0030s ; Time = 1.97s ; ETA = 5.04s Iter = 808 ; Dt = 0.0003s ; Time = 2.01s ; ETA = 0.00s =================================================================== | END | =================================================================== Problem.solve() completed in 307.40450263023376s (wall time) .. GENERATED FROM PYTHON SOURCE LINES 262-265 Errors plots: -------------------- .. GENERATED FROM PYTHON SOURCE LINES 266-315 .. code-block:: default def plot_error(ax, error): err = np.zeros((N_mesh)) ab = np.zeros((N_mesh)) if SECOND_ORDER: ORD = 2 else: ORD = 1 for order in range(ORD): for M in range(N_mesh): err[M] = np.log(error[M, order]) #ab[M] = 0.5*np.log(NT[M]/NT[0]) ab[M] = 0.5*np.log((NT[M]*NL[M])/(NT[0]*NL[0])) ax.plot(ab, err, color=colors[order], marker=markers[order], label=labels[order]) def plot_reference(ax, error): xmax = 0.5*np.log(NT[N_mesh-1]/NT[0]) ax.plot([0.0, xmax], [np.log(error[0, 0]), -1.0*xmax+np.log(error[0, 0])],\ color='k', ls=':', lw=0.5) ax.plot([0.0, xmax], [np.log(error[0, 1]), -2.0*xmax+np.log(error[0, 1])],\ color='k', ls='-.', lw=0.5) if PLOT_ERROR == True: plt.rcParams["figure.figsize"] = [10, 5] colors = ['blue','red','green','purple'] markers = ['s','*','o','v'] labels = ['$1^{st}$ order','$2^{nd}$ order'] fig = plt.figure() ax3 = fig.add_subplot(111) plot_error(ax3, error_Rho_cumul) plot_reference(ax3, error_Rho_cumul) ax3.set_xlabel('$log(l_0/l_i)$') ax3.set_ylabel('Error of $\\rho h$ in $L_2$ norm') ax3.set_title("Convergence of $\\rho h$") plt.legend(loc=1) plt.show() fig = plt.figure() ax0 = fig.add_subplot(111) plot_error(ax0, error_H_cumul) plot_reference(ax0, error_H_cumul) ax0.set_xlabel('$log(l_0/l_i)$') ax0.set_ylabel('Error of $h$ in $L_2$ norm') ax0.set_title("Convergence of $h$") plt.legend(loc=1) plt.show() .. rst-class:: sphx-glr-horizontal * .. image-sg:: /auto_examples_mesh/images/sphx_glr_example_meshconvergence_thacker3d_006.png :alt: Convergence of $\rho h$ :srcset: /auto_examples_mesh/images/sphx_glr_example_meshconvergence_thacker3d_006.png :class: sphx-glr-multi-img * .. image-sg:: /auto_examples_mesh/images/sphx_glr_example_meshconvergence_thacker3d_007.png :alt: Convergence of $h$ :srcset: /auto_examples_mesh/images/sphx_glr_example_meshconvergence_thacker3d_007.png :class: sphx-glr-multi-img .. rst-class:: sphx-glr-timing **Total running time of the script:** ( 11 minutes 4.242 seconds) .. _sphx_glr_download_auto_examples_mesh_example_meshconvergence_thacker3d.py: .. only:: html .. container:: sphx-glr-footer sphx-glr-footer-example .. container:: sphx-glr-download sphx-glr-download-python :download:`Download Python source code: example_meshconvergence_thacker3d.py ` .. container:: sphx-glr-download sphx-glr-download-jupyter :download:`Download Jupyter notebook: example_meshconvergence_thacker3d.ipynb ` .. only:: html .. rst-class:: sphx-glr-signature `Gallery generated by Sphinx-Gallery `_