.. DO NOT EDIT. .. THIS FILE WAS AUTOMATICALLY GENERATED BY SPHINX-GALLERY. .. TO MAKE CHANGES, EDIT THE SOURCE PYTHON FILE: .. "auto_examples_mesh/example_meshconvergence_bump.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_bump.py: ================================================================================ Bump mesh convergence ================================================================================ In this example mesh convergence is carried out in the case of a bump. H and QX errors are computed by comparison to single layer 1D analytical solution. .. GENERATED FROM PYTHON SOURCE LINES 10-21 .. code-block:: default import os import h5py import numpy as np import matplotlib.pyplot as plt import freshkiss3d as fk import freshkiss3d.extra.plots as fk_plt #sphinx_gallery_thumbnail_number = 5 os.system('rm -r outputs') .. rst-class:: sphx-glr-script-out .. code-block:: none 256 .. GENERATED FROM PYTHON SOURCE LINES 22-24 Mesh parameters: -------------------- .. GENERATED FROM PYTHON SOURCE LINES 25-30 .. code-block:: default N_mesh = 4 #Number of mesh for mesh convergence First_mesh = 2 NL = [1 for r in range(N_mesh)] .. GENERATED FROM PYTHON SOURCE LINES 31-33 Case parameters: -------------------- .. GENERATED FROM PYTHON SOURCE LINES 34-69 .. code-block:: default FLOW = 'sub' if FLOW == 'sub': FREE_SURFACE_0 = 0.7 Q_IN = 0.2 if FLOW == 'trans': FREE_SURFACE_0 = 0.7 Q_IN = 1.5 if FLOW == 'trans_shock': FREE_SURFACE_0 = 0.33 Q_IN = 0.18 H_OUT = FREE_SURFACE_0 FINAL_TIME = 1000. WHEN = [5.99*FINAL_TIME/6.] MESH_SPLITTING = False PLOT_MESH = True PLOT_SOL = False PLOT_ERROR = True SAVE_ERROR = True X_B = 10. def topo(x): if 8.< x <12.: topo = 0.2 - 0.05*(x - X_B)**2 else: topo = 0.0 return topo def topo_gaussian(x): topo = 0.25*np.exp(-0.5*(x-X_B)**2) return topo .. GENERATED FROM PYTHON SOURCE LINES 70-72 Lists initialization: ----------------------- .. GENERATED FROM PYTHON SOURCE LINES 73-83 .. code-block:: default cases = [] simutime_list = [] triangular_mesh_list = [] NT = [] layer_list = [] primitives_list = [] problem_list = [] bump_analytic_list = [] .. GENERATED FROM PYTHON SOURCE LINES 84-86 Meshes: -------------------- .. GENERATED FROM PYTHON SOURCE LINES 87-110 .. code-block:: default for M in range(N_mesh): if MESH_SPLITTING: #Generates first mesh from bump0.geo and the other by splitting of the first mesh if M == 0: os.system('gmsh -2 -format msh2 inputs/bump{}.geo -o inputs/bump{}.msh'.format(First_mesh, First_mesh)) TG, vertex_labels, boundary_labels = fk.read_msh('inputs/bump{}.msh'.format(First_mesh)) triangular_mesh = fk.TriangularMesh(TG, vertex_labels, boundary_labels) triangular_mesh_list.append(triangular_mesh) os.system('rm inputs/bump{}.msh'.format(First_mesh)) else: triangular_mesh = triangular_mesh.refine_by_splitting() triangular_mesh_list.append(triangular_mesh) else: #Generates all meshes from bump#.geo files os.system('gmsh -2 -format msh2 inputs/bump{}.geo -o inputs/bump{}.msh'.format(M, M)) TG, vertex_labels, boundary_labels = fk.read_msh('inputs/bump{}.msh'.format(M)) triangular_mesh_list.append(fk.TriangularMesh(TG, vertex_labels, boundary_labels)) os.system('rm inputs/bump{}.msh'.format(M)) 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_bump_001.png :alt: Mesh :srcset: /auto_examples_mesh/images/sphx_glr_example_meshconvergence_bump_001.png :class: sphx-glr-multi-img * .. image-sg:: /auto_examples_mesh/images/sphx_glr_example_meshconvergence_bump_002.png :alt: Mesh :srcset: /auto_examples_mesh/images/sphx_glr_example_meshconvergence_bump_002.png :class: sphx-glr-multi-img * .. image-sg:: /auto_examples_mesh/images/sphx_glr_example_meshconvergence_bump_003.png :alt: Mesh :srcset: /auto_examples_mesh/images/sphx_glr_example_meshconvergence_bump_003.png :class: sphx-glr-multi-img * .. image-sg:: /auto_examples_mesh/images/sphx_glr_example_meshconvergence_bump_004.png :alt: Mesh :srcset: /auto_examples_mesh/images/sphx_glr_example_meshconvergence_bump_004.png :class: sphx-glr-multi-img .. GENERATED FROM PYTHON SOURCE LINES 111-116 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 117-134 .. 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 = False if space_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)) cases.append(case(NL=NL[M], mesh_id=M, time_order=True, space_order=True)) .. GENERATED FROM PYTHON SOURCE LINES 135-139 Boundary conditions -------------------- Here we define boundary conditions, which are the same for every cases. .. GENERATED FROM PYTHON SOURCE LINES 140-148 .. code-block:: default fluvial_flowrates = [fk.FluvialFlowRate(ref=4, flowrate=Q_IN, x_flux_direction=1.0, y_flux_direction=0.0)] fluvial_heights = [fk.FluvialHeight(ref=2, height=H_OUT)] torrential_outs = [fk.TorrentialOut(ref=2)] slides = [fk.Slide(ref=3), fk.Slide(ref=1)] .. GENERATED FROM PYTHON SOURCE LINES 149-151 Solving (loop over all cases): ------------------------------- .. GENERATED FROM PYTHON SOURCE LINES 152-219 .. code-block:: default for I, case in enumerate(cases): print(" ") print(" ") print(" SOLVING CASE {} ".format(I)) print(" ") print(" ") mesh_id = case.mesh_id time_second_order = case.time_order numerical_params = { 'space_second_order':case.space_order, 'ipres':case.ipres} simutime_list.append(fk.SimuTime( final_time=FINAL_TIME, time_iteration_max=1000000, second_order=time_second_order)) bump_analytic_list.append(fk.Bump( triangular_mesh_list[mesh_id], case=FLOW, q_in=Q_IN, h_out=H_OUT, free_surface_init=FREE_SURFACE_0, x_b = X_B, topo = topo_gaussian, scheduler=fk.schedules(times=WHEN), compute_error=True, error_type='L2', error_output='txt')) bump_analytic_list[I](0.) layer_list.append(fk.Layer( case.NL, triangular_mesh_list[mesh_id], topography=bump_analytic_list[I].topography)) primitives_list.append(fk.Primitive( triangular_mesh_list[mesh_id], layer_list[I], free_surface=FREE_SURFACE_0, QXinit=Q_IN, QYinit=0.)) if PLOT_SOL: create_figure = {'plot':fk_plt.plot_freesurface_3d_analytic} create_figure_scheduler = fk.schedules(times=[0.99*FINAL_TIME]) else: create_figure = None create_figure_scheduler = None problem_list.append(fk.Problem( simutime_list[I], triangular_mesh_list[mesh_id], layer_list[I], primitives_list[I], analytic_sol=bump_analytic_list[I], numerical_parameters=numerical_params, fluvial_flowrates=fluvial_flowrates, fluvial_heights=fluvial_heights, slides=slides, custom_funct=create_figure, custom_funct_scheduler=create_figure_scheduler)) problem_list[I].solve() if PLOT_SOL: plt.show() .. rst-class:: sphx-glr-script-out .. code-block:: none SOLVING CASE 0 =================================================================== | INITIALIZATION | =================================================================== Problem size: Nodes=74, Layers=1, Triangles=96, Iter = 0 ; Dt = 0.0000s ; Time = 0.00s ; ETA = 0.00s =================================================================== | TIME LOOP | =================================================================== Iter = 724 ; Dt = 0.0282s ; Time = 20.43s ; ETA = 6.65s Iter = 1446 ; Dt = 0.0283s ; Time = 40.82s ; ETA = 6.41s Iter = 2169 ; Dt = 0.0282s ; Time = 61.23s ; ETA = 6.54s Iter = 2892 ; Dt = 0.0282s ; Time = 81.64s ; ETA = 6.19s Iter = 3615 ; Dt = 0.0282s ; Time = 102.05s ; ETA = 6.47s Iter = 4338 ; Dt = 0.0282s ; Time = 122.46s ; ETA = 5.85s Iter = 5061 ; Dt = 0.0282s ; Time = 142.87s ; ETA = 5.76s Iter = 5784 ; Dt = 0.0282s ; Time = 163.28s ; ETA = 5.69s Iter = 6507 ; Dt = 0.0282s ; Time = 183.69s ; ETA = 5.60s Iter = 7230 ; Dt = 0.0282s ; Time = 204.10s ; ETA = 5.33s Iter = 7953 ; Dt = 0.0282s ; Time = 224.51s ; ETA = 5.23s Iter = 8676 ; Dt = 0.0282s ; Time = 244.92s ; ETA = 5.36s Iter = 9398 ; Dt = 0.0282s ; Time = 265.31s ; ETA = 4.91s Iter = 10121 ; Dt = 0.0282s ; Time = 285.72s ; ETA = 4.81s Iter = 10844 ; Dt = 0.0282s ; Time = 306.13s ; ETA = 4.86s Iter = 11567 ; Dt = 0.0282s ; Time = 326.54s ; ETA = 4.57s Iter = 12290 ; Dt = 0.0282s ; Time = 346.95s ; ETA = 4.65s Iter = 13013 ; Dt = 0.0282s ; Time = 367.36s ; ETA = 4.48s Iter = 13736 ; Dt = 0.0282s ; Time = 387.77s ; ETA = 4.28s Iter = 14459 ; Dt = 0.0282s ; Time = 408.18s ; ETA = 3.94s Iter = 15182 ; Dt = 0.0282s ; Time = 428.59s ; ETA = 3.83s Iter = 15905 ; Dt = 0.0282s ; Time = 449.00s ; ETA = 4.03s Iter = 16628 ; Dt = 0.0282s ; Time = 469.41s ; ETA = 3.58s Iter = 17350 ; Dt = 0.0282s ; Time = 489.80s ; ETA = 4.10s Iter = 18073 ; Dt = 0.0282s ; Time = 510.21s ; ETA = 3.38s Iter = 18796 ; Dt = 0.0282s ; Time = 530.62s ; ETA = 3.38s Iter = 19519 ; Dt = 0.0282s ; Time = 551.03s ; ETA = 3.23s Iter = 20242 ; Dt = 0.0282s ; Time = 571.44s ; ETA = 3.08s Iter = 20965 ; Dt = 0.0282s ; Time = 591.85s ; ETA = 2.73s Iter = 21688 ; Dt = 0.0282s ; Time = 612.26s ; ETA = 2.62s Iter = 22411 ; Dt = 0.0282s ; Time = 632.67s ; ETA = 2.50s Iter = 23134 ; Dt = 0.0282s ; Time = 653.08s ; ETA = 2.35s Iter = 23857 ; Dt = 0.0282s ; Time = 673.49s ; ETA = 2.19s Iter = 24580 ; Dt = 0.0282s ; Time = 693.90s ; ETA = 2.08s Iter = 25302 ; Dt = 0.0282s ; Time = 714.29s ; ETA = 1.92s Iter = 26025 ; Dt = 0.0282s ; Time = 734.70s ; ETA = 1.91s Iter = 26748 ; Dt = 0.0282s ; Time = 755.11s ; ETA = 1.64s Iter = 27471 ; Dt = 0.0282s ; Time = 775.52s ; ETA = 1.50s Iter = 28194 ; Dt = 0.0282s ; Time = 795.93s ; ETA = 1.36s Iter = 28917 ; Dt = 0.0282s ; Time = 816.34s ; ETA = 1.24s Iter = 29640 ; Dt = 0.0282s ; Time = 836.75s ; ETA = 1.11s Iter = 30363 ; Dt = 0.0282s ; Time = 857.16s ; ETA = 1.01s Iter = 31086 ; Dt = 0.0282s ; Time = 877.57s ; ETA = 0.82s Iter = 31809 ; Dt = 0.0282s ; Time = 897.98s ; ETA = 0.68s Iter = 32532 ; Dt = 0.0282s ; Time = 918.40s ; ETA = 0.55s Iter = 33254 ; Dt = 0.0282s ; Time = 938.78s ; ETA = 0.41s Iter = 33977 ; Dt = 0.0282s ; Time = 959.19s ; ETA = 0.28s Iter = 34700 ; Dt = 0.0282s ; Time = 979.60s ; ETA = 0.15s Iter = 35423 ; Dt = 0.0182s ; Time = 1000.00s ; ETA = 0.00s =================================================================== | END | =================================================================== Problem.solve() completed in 6.954565525054932s (wall time) SOLVING CASE 1 =================================================================== | INITIALIZATION | =================================================================== Problem size: Nodes=74, Layers=1, Triangles=96, Iter = 0 ; Dt = 0.0000s ; Time = 0.00s ; ETA = 0.00s =================================================================== | TIME LOOP | =================================================================== Iter = 722 ; Dt = 0.0283s ; Time = 20.41s ; ETA = 16.62s Iter = 1444 ; Dt = 0.0283s ; Time = 40.83s ; ETA = 16.19s Iter = 2166 ; Dt = 0.0283s ; Time = 61.25s ; ETA = 15.85s Iter = 2887 ; Dt = 0.0283s ; Time = 81.64s ; ETA = 15.47s Iter = 3609 ; Dt = 0.0283s ; Time = 102.06s ; ETA = 15.18s Iter = 4330 ; Dt = 0.0283s ; Time = 122.45s ; ETA = 14.82s Iter = 5052 ; Dt = 0.0283s ; Time = 142.87s ; ETA = 14.37s Iter = 5774 ; Dt = 0.0283s ; Time = 163.29s ; ETA = 14.36s Iter = 6495 ; Dt = 0.0283s ; Time = 183.68s ; ETA = 13.38s Iter = 7217 ; Dt = 0.0283s ; Time = 204.11s ; ETA = 13.52s Iter = 7938 ; Dt = 0.0283s ; Time = 224.50s ; ETA = 12.87s Iter = 8660 ; Dt = 0.0283s ; Time = 244.92s ; ETA = 12.41s Iter = 9381 ; Dt = 0.0283s ; Time = 265.31s ; ETA = 12.10s Iter = 10103 ; Dt = 0.0283s ; Time = 285.73s ; ETA = 12.14s Iter = 10824 ; Dt = 0.0283s ; Time = 306.12s ; ETA = 11.61s Iter = 11546 ; Dt = 0.0283s ; Time = 326.54s ; ETA = 11.92s Iter = 12268 ; Dt = 0.0283s ; Time = 346.96s ; ETA = 10.68s Iter = 12989 ; Dt = 0.0283s ; Time = 367.36s ; ETA = 10.70s Iter = 13711 ; Dt = 0.0283s ; Time = 387.78s ; ETA = 10.03s Iter = 14432 ; Dt = 0.0283s ; Time = 408.17s ; ETA = 9.72s Iter = 15154 ; Dt = 0.0283s ; Time = 428.59s ; ETA = 9.39s Iter = 15875 ; Dt = 0.0283s ; Time = 448.98s ; ETA = 9.41s Iter = 16597 ; Dt = 0.0283s ; Time = 469.40s ; ETA = 8.96s Iter = 17319 ; Dt = 0.0283s ; Time = 489.82s ; ETA = 8.70s Iter = 18040 ; Dt = 0.0283s ; Time = 510.21s ; ETA = 8.27s Iter = 18762 ; Dt = 0.0283s ; Time = 530.64s ; ETA = 7.69s Iter = 19483 ; Dt = 0.0283s ; Time = 551.03s ; ETA = 7.52s Iter = 20205 ; Dt = 0.0283s ; Time = 571.45s ; ETA = 7.26s Iter = 20926 ; Dt = 0.0283s ; Time = 591.84s ; ETA = 6.71s Iter = 21648 ; Dt = 0.0283s ; Time = 612.26s ; ETA = 6.38s Iter = 22370 ; Dt = 0.0283s ; Time = 632.68s ; ETA = 6.04s Iter = 23091 ; Dt = 0.0283s ; Time = 653.07s ; ETA = 5.73s Iter = 23813 ; Dt = 0.0283s ; Time = 673.49s ; ETA = 5.58s Iter = 24534 ; Dt = 0.0283s ; Time = 693.89s ; ETA = 5.18s Iter = 25256 ; Dt = 0.0283s ; Time = 714.31s ; ETA = 4.88s Iter = 25977 ; Dt = 0.0283s ; Time = 734.70s ; ETA = 4.35s Iter = 26699 ; Dt = 0.0283s ; Time = 755.12s ; ETA = 4.18s Iter = 27420 ; Dt = 0.0283s ; Time = 775.51s ; ETA = 3.88s Iter = 28142 ; Dt = 0.0283s ; Time = 795.93s ; ETA = 3.31s Iter = 28864 ; Dt = 0.0283s ; Time = 816.35s ; ETA = 3.26s Iter = 29585 ; Dt = 0.0283s ; Time = 836.74s ; ETA = 2.77s Iter = 30307 ; Dt = 0.0283s ; Time = 857.16s ; ETA = 2.39s Iter = 31028 ; Dt = 0.0283s ; Time = 877.56s ; ETA = 2.08s Iter = 31750 ; Dt = 0.0283s ; Time = 897.98s ; ETA = 1.65s Iter = 32471 ; Dt = 0.0283s ; Time = 918.37s ; ETA = 1.39s Iter = 33193 ; Dt = 0.0283s ; Time = 938.79s ; ETA = 1.03s Iter = 33915 ; Dt = 0.0283s ; Time = 959.21s ; ETA = 0.67s Iter = 34636 ; Dt = 0.0283s ; Time = 979.60s ; ETA = 0.36s Iter = 35358 ; Dt = 0.0050s ; Time = 1000.00s ; ETA = 0.00s =================================================================== | END | =================================================================== Problem.solve() completed in 16.951236248016357s (wall time) SOLVING CASE 2 =================================================================== | INITIALIZATION | =================================================================== Problem size: Nodes=149, Layers=1, Triangles=196, Iter = 0 ; Dt = 0.0000s ; Time = 0.00s ; ETA = 0.00s =================================================================== | TIME LOOP | =================================================================== Iter = 997 ; Dt = 0.0204s ; Time = 20.41s ; ETA = 10.37s Iter = 1994 ; Dt = 0.0205s ; Time = 40.83s ; ETA = 10.01s Iter = 2990 ; Dt = 0.0205s ; Time = 61.23s ; ETA = 10.42s Iter = 3987 ; Dt = 0.0205s ; Time = 81.65s ; ETA = 10.74s Iter = 4983 ; Dt = 0.0205s ; Time = 102.04s ; ETA = 9.68s Iter = 5980 ; Dt = 0.0205s ; Time = 122.46s ; ETA = 9.15s Iter = 6977 ; Dt = 0.0205s ; Time = 142.88s ; ETA = 10.48s Iter = 7973 ; Dt = 0.0205s ; Time = 163.28s ; ETA = 8.72s Iter = 8970 ; Dt = 0.0205s ; Time = 183.69s ; ETA = 8.58s Iter = 9966 ; Dt = 0.0205s ; Time = 204.09s ; ETA = 8.50s Iter = 10963 ; Dt = 0.0205s ; Time = 224.51s ; ETA = 8.08s Iter = 11959 ; Dt = 0.0205s ; Time = 244.90s ; ETA = 7.92s Iter = 12956 ; Dt = 0.0205s ; Time = 265.32s ; ETA = 7.65s Iter = 13952 ; Dt = 0.0205s ; Time = 285.72s ; ETA = 7.49s Iter = 14949 ; Dt = 0.0205s ; Time = 306.14s ; ETA = 7.67s Iter = 15945 ; Dt = 0.0205s ; Time = 326.53s ; ETA = 7.13s Iter = 16942 ; Dt = 0.0205s ; Time = 346.95s ; ETA = 6.84s Iter = 17938 ; Dt = 0.0205s ; Time = 367.35s ; ETA = 6.61s Iter = 18935 ; Dt = 0.0205s ; Time = 387.77s ; ETA = 7.04s Iter = 19931 ; Dt = 0.0205s ; Time = 408.16s ; ETA = 6.56s Iter = 20928 ; Dt = 0.0205s ; Time = 428.58s ; ETA = 7.22s Iter = 21925 ; Dt = 0.0205s ; Time = 449.00s ; ETA = 6.14s Iter = 22921 ; Dt = 0.0205s ; Time = 469.40s ; ETA = 5.58s Iter = 23918 ; Dt = 0.0205s ; Time = 489.81s ; ETA = 5.64s Iter = 24914 ; Dt = 0.0205s ; Time = 510.21s ; ETA = 5.10s Iter = 25911 ; Dt = 0.0205s ; Time = 530.63s ; ETA = 4.92s Iter = 26907 ; Dt = 0.0205s ; Time = 551.03s ; ETA = 4.74s Iter = 27904 ; Dt = 0.0205s ; Time = 571.44s ; ETA = 4.52s Iter = 28900 ; Dt = 0.0205s ; Time = 591.84s ; ETA = 4.30s Iter = 29897 ; Dt = 0.0205s ; Time = 612.26s ; ETA = 4.29s Iter = 30893 ; Dt = 0.0205s ; Time = 632.66s ; ETA = 4.38s Iter = 31890 ; Dt = 0.0205s ; Time = 653.07s ; ETA = 3.70s Iter = 32886 ; Dt = 0.0205s ; Time = 673.47s ; ETA = 3.60s Iter = 33883 ; Dt = 0.0205s ; Time = 693.89s ; ETA = 3.19s Iter = 34880 ; Dt = 0.0205s ; Time = 714.31s ; ETA = 3.09s Iter = 35876 ; Dt = 0.0205s ; Time = 734.70s ; ETA = 2.77s Iter = 36873 ; Dt = 0.0205s ; Time = 755.12s ; ETA = 2.57s Iter = 37869 ; Dt = 0.0205s ; Time = 775.52s ; ETA = 2.35s Iter = 38866 ; Dt = 0.0205s ; Time = 795.94s ; ETA = 2.39s Iter = 39862 ; Dt = 0.0205s ; Time = 816.33s ; ETA = 1.93s Iter = 40859 ; Dt = 0.0205s ; Time = 836.75s ; ETA = 1.73s Iter = 41855 ; Dt = 0.0205s ; Time = 857.15s ; ETA = 1.50s Iter = 42852 ; Dt = 0.0205s ; Time = 877.57s ; ETA = 1.28s Iter = 43848 ; Dt = 0.0205s ; Time = 897.96s ; ETA = 1.12s Iter = 44845 ; Dt = 0.0205s ; Time = 918.38s ; ETA = 0.85s Iter = 45841 ; Dt = 0.0205s ; Time = 938.78s ; ETA = 0.64s Iter = 46838 ; Dt = 0.0205s ; Time = 959.19s ; ETA = 0.43s Iter = 47834 ; Dt = 0.0205s ; Time = 979.59s ; ETA = 0.22s Iter = 48831 ; Dt = 0.0109s ; Time = 1000.00s ; ETA = 0.00s =================================================================== | END | =================================================================== Problem.solve() completed in 10.847069263458252s (wall time) SOLVING CASE 3 =================================================================== | INITIALIZATION | =================================================================== Problem size: Nodes=149, Layers=1, Triangles=196, Iter = 0 ; Dt = 0.0000s ; Time = 0.00s ; ETA = 0.00s =================================================================== | TIME LOOP | =================================================================== Iter = 996 ; Dt = 0.0205s ; Time = 20.42s ; ETA = 32.17s Iter = 1991 ; Dt = 0.0205s ; Time = 40.84s ; ETA = 32.20s Iter = 2986 ; Dt = 0.0205s ; Time = 61.24s ; ETA = 31.38s Iter = 3980 ; Dt = 0.0205s ; Time = 81.63s ; ETA = 29.50s Iter = 4975 ; Dt = 0.0205s ; Time = 102.04s ; ETA = 30.04s Iter = 5970 ; Dt = 0.0205s ; Time = 122.46s ; ETA = 30.07s Iter = 6965 ; Dt = 0.0205s ; Time = 142.87s ; ETA = 30.40s Iter = 7960 ; Dt = 0.0205s ; Time = 163.28s ; ETA = 27.05s Iter = 8955 ; Dt = 0.0205s ; Time = 183.69s ; ETA = 28.15s Iter = 9949 ; Dt = 0.0205s ; Time = 204.08s ; ETA = 26.19s Iter = 10944 ; Dt = 0.0205s ; Time = 224.50s ; ETA = 25.84s Iter = 11939 ; Dt = 0.0205s ; Time = 244.91s ; ETA = 24.87s Iter = 12934 ; Dt = 0.0205s ; Time = 265.32s ; ETA = 24.70s Iter = 13929 ; Dt = 0.0205s ; Time = 285.73s ; ETA = 23.55s Iter = 14923 ; Dt = 0.0205s ; Time = 306.12s ; ETA = 23.40s Iter = 15918 ; Dt = 0.0205s ; Time = 326.54s ; ETA = 22.37s Iter = 16913 ; Dt = 0.0205s ; Time = 346.95s ; ETA = 24.17s Iter = 17908 ; Dt = 0.0205s ; Time = 367.36s ; ETA = 20.53s Iter = 18903 ; Dt = 0.0205s ; Time = 387.77s ; ETA = 20.29s Iter = 19898 ; Dt = 0.0205s ; Time = 408.18s ; ETA = 19.55s Iter = 20892 ; Dt = 0.0205s ; Time = 428.57s ; ETA = 18.87s Iter = 21887 ; Dt = 0.0205s ; Time = 448.99s ; ETA = 18.19s Iter = 22882 ; Dt = 0.0205s ; Time = 469.40s ; ETA = 17.92s Iter = 23877 ; Dt = 0.0205s ; Time = 489.81s ; ETA = 17.14s Iter = 24872 ; Dt = 0.0205s ; Time = 510.22s ; ETA = 18.19s Iter = 25866 ; Dt = 0.0205s ; Time = 530.61s ; ETA = 16.00s Iter = 26861 ; Dt = 0.0205s ; Time = 551.03s ; ETA = 16.48s Iter = 27856 ; Dt = 0.0205s ; Time = 571.44s ; ETA = 15.20s Iter = 28851 ; Dt = 0.0205s ; Time = 591.85s ; ETA = 13.43s Iter = 29846 ; Dt = 0.0205s ; Time = 612.26s ; ETA = 12.67s Iter = 30840 ; Dt = 0.0205s ; Time = 632.65s ; ETA = 12.13s Iter = 31835 ; Dt = 0.0205s ; Time = 653.07s ; ETA = 11.33s Iter = 32830 ; Dt = 0.0205s ; Time = 673.48s ; ETA = 10.77s Iter = 33825 ; Dt = 0.0205s ; Time = 693.89s ; ETA = 10.46s Iter = 34820 ; Dt = 0.0205s ; Time = 714.30s ; ETA = 9.19s Iter = 35815 ; Dt = 0.0205s ; Time = 734.71s ; ETA = 8.95s Iter = 36809 ; Dt = 0.0205s ; Time = 755.10s ; ETA = 8.01s Iter = 37804 ; Dt = 0.0205s ; Time = 775.52s ; ETA = 7.26s Iter = 38799 ; Dt = 0.0205s ; Time = 795.93s ; ETA = 6.75s Iter = 39794 ; Dt = 0.0205s ; Time = 816.34s ; ETA = 6.15s Iter = 40789 ; Dt = 0.0205s ; Time = 836.75s ; ETA = 5.41s Iter = 41783 ; Dt = 0.0205s ; Time = 857.14s ; ETA = 4.61s Iter = 42778 ; Dt = 0.0205s ; Time = 877.56s ; ETA = 4.10s Iter = 43773 ; Dt = 0.0205s ; Time = 897.97s ; ETA = 3.28s Iter = 44768 ; Dt = 0.0205s ; Time = 918.38s ; ETA = 2.69s Iter = 45763 ; Dt = 0.0205s ; Time = 938.79s ; ETA = 2.09s Iter = 46758 ; Dt = 0.0205s ; Time = 959.20s ; ETA = 1.41s Iter = 47752 ; Dt = 0.0205s ; Time = 979.60s ; ETA = 1.39s Iter = 48747 ; Dt = 0.0132s ; Time = 1000.00s ; ETA = 0.00s =================================================================== | END | =================================================================== Problem.solve() completed in 33.620768308639526s (wall time) SOLVING CASE 4 =================================================================== | INITIALIZATION | =================================================================== Problem size: Nodes=346, Layers=1, Triangles=548, Iter = 0 ; Dt = 0.0000s ; Time = 0.00s ; ETA = 0.00s =================================================================== | TIME LOOP | =================================================================== Iter = 1647 ; Dt = 0.0124s ; Time = 20.41s ; ETA = 25.55s Iter = 3293 ; Dt = 0.0124s ; Time = 40.82s ; ETA = 23.08s Iter = 4939 ; Dt = 0.0124s ; Time = 61.23s ; ETA = 24.22s Iter = 6585 ; Dt = 0.0124s ; Time = 81.63s ; ETA = 23.20s Iter = 8232 ; Dt = 0.0124s ; Time = 102.05s ; ETA = 22.88s Iter = 9878 ; Dt = 0.0124s ; Time = 122.46s ; ETA = 21.30s Iter = 11524 ; Dt = 0.0124s ; Time = 142.87s ; ETA = 21.40s Iter = 13170 ; Dt = 0.0124s ; Time = 163.27s ; ETA = 20.27s Iter = 14816 ; Dt = 0.0124s ; Time = 183.68s ; ETA = 20.53s Iter = 16462 ; Dt = 0.0124s ; Time = 204.09s ; ETA = 19.48s Iter = 18108 ; Dt = 0.0124s ; Time = 224.49s ; ETA = 19.60s Iter = 19754 ; Dt = 0.0124s ; Time = 244.90s ; ETA = 18.56s Iter = 21400 ; Dt = 0.0124s ; Time = 265.31s ; ETA = 17.75s Iter = 23047 ; Dt = 0.0124s ; Time = 285.73s ; ETA = 18.02s Iter = 24693 ; Dt = 0.0124s ; Time = 306.13s ; ETA = 17.21s Iter = 26339 ; Dt = 0.0124s ; Time = 326.54s ; ETA = 16.40s Iter = 27985 ; Dt = 0.0124s ; Time = 346.95s ; ETA = 16.40s Iter = 29631 ; Dt = 0.0124s ; Time = 367.35s ; ETA = 15.35s Iter = 31277 ; Dt = 0.0124s ; Time = 387.76s ; ETA = 15.35s Iter = 32923 ; Dt = 0.0124s ; Time = 408.17s ; ETA = 14.31s Iter = 34569 ; Dt = 0.0124s ; Time = 428.57s ; ETA = 14.49s Iter = 36216 ; Dt = 0.0124s ; Time = 448.99s ; ETA = 14.29s Iter = 37862 ; Dt = 0.0124s ; Time = 469.40s ; ETA = 13.48s Iter = 39508 ; Dt = 0.0124s ; Time = 489.81s ; ETA = 12.37s Iter = 41154 ; Dt = 0.0124s ; Time = 510.21s ; ETA = 12.53s Iter = 42800 ; Dt = 0.0124s ; Time = 530.62s ; ETA = 12.56s Iter = 44446 ; Dt = 0.0124s ; Time = 551.03s ; ETA = 11.01s Iter = 46092 ; Dt = 0.0124s ; Time = 571.43s ; ETA = 10.39s Iter = 47738 ; Dt = 0.0124s ; Time = 591.84s ; ETA = 9.94s Iter = 49384 ; Dt = 0.0124s ; Time = 612.25s ; ETA = 9.43s Iter = 51031 ; Dt = 0.0124s ; Time = 632.66s ; ETA = 9.36s Iter = 52677 ; Dt = 0.0124s ; Time = 653.07s ; ETA = 8.44s Iter = 54323 ; Dt = 0.0124s ; Time = 673.48s ; ETA = 8.06s Iter = 55969 ; Dt = 0.0124s ; Time = 693.88s ; ETA = 7.42s Iter = 57615 ; Dt = 0.0124s ; Time = 714.29s ; ETA = 6.93s Iter = 59261 ; Dt = 0.0124s ; Time = 734.70s ; ETA = 6.43s Iter = 60907 ; Dt = 0.0124s ; Time = 755.11s ; ETA = 6.54s Iter = 62553 ; Dt = 0.0124s ; Time = 775.51s ; ETA = 5.43s Iter = 64199 ; Dt = 0.0124s ; Time = 795.92s ; ETA = 4.92s Iter = 65846 ; Dt = 0.0124s ; Time = 816.34s ; ETA = 4.80s Iter = 67492 ; Dt = 0.0124s ; Time = 836.74s ; ETA = 3.95s Iter = 69138 ; Dt = 0.0124s ; Time = 857.15s ; ETA = 3.60s Iter = 70784 ; Dt = 0.0124s ; Time = 877.56s ; ETA = 2.98s Iter = 72430 ; Dt = 0.0124s ; Time = 897.96s ; ETA = 2.57s Iter = 74076 ; Dt = 0.0124s ; Time = 918.37s ; ETA = 1.98s Iter = 75722 ; Dt = 0.0124s ; Time = 938.78s ; ETA = 1.55s Iter = 77368 ; Dt = 0.0124s ; Time = 959.18s ; ETA = 0.99s Iter = 79015 ; Dt = 0.0124s ; Time = 979.60s ; ETA = 0.51s Iter = 80661 ; Dt = 0.0018s ; Time = 1000.00s ; ETA = 0.00s =================================================================== | END | =================================================================== Problem.solve() completed in 25.04877996444702s (wall time) SOLVING CASE 5 =================================================================== | INITIALIZATION | =================================================================== Problem size: Nodes=346, Layers=1, Triangles=548, Iter = 0 ; Dt = 0.0000s ; Time = 0.00s ; ETA = 0.00s =================================================================== | TIME LOOP | =================================================================== Iter = 1733 ; Dt = 0.0118s ; Time = 20.41s ; ETA = 103.55s Iter = 3465 ; Dt = 0.0118s ; Time = 40.82s ; ETA = 103.87s Iter = 5198 ; Dt = 0.0118s ; Time = 61.23s ; ETA = 100.92s Iter = 6930 ; Dt = 0.0118s ; Time = 81.64s ; ETA = 99.60s Iter = 8663 ; Dt = 0.0118s ; Time = 102.05s ; ETA = 104.54s Iter = 10395 ; Dt = 0.0118s ; Time = 122.46s ; ETA = 97.47s Iter = 12127 ; Dt = 0.0118s ; Time = 142.86s ; ETA = 92.65s Iter = 13859 ; Dt = 0.0118s ; Time = 163.27s ; ETA = 91.67s Iter = 15591 ; Dt = 0.0118s ; Time = 183.68s ; ETA = 88.81s Iter = 17323 ; Dt = 0.0118s ; Time = 204.08s ; ETA = 84.27s Iter = 19056 ; Dt = 0.0118s ; Time = 224.50s ; ETA = 83.80s Iter = 20788 ; Dt = 0.0118s ; Time = 244.91s ; ETA = 79.87s Iter = 22520 ; Dt = 0.0118s ; Time = 265.31s ; ETA = 80.53s Iter = 24252 ; Dt = 0.0118s ; Time = 285.72s ; ETA = 77.32s Iter = 25984 ; Dt = 0.0118s ; Time = 306.13s ; ETA = 74.70s Iter = 27716 ; Dt = 0.0118s ; Time = 326.53s ; ETA = 72.50s Iter = 29448 ; Dt = 0.0118s ; Time = 346.94s ; ETA = 69.70s Iter = 31181 ; Dt = 0.0118s ; Time = 367.36s ; ETA = 67.15s Iter = 32913 ; Dt = 0.0118s ; Time = 387.76s ; ETA = 73.85s Iter = 34645 ; Dt = 0.0118s ; Time = 408.17s ; ETA = 71.77s Iter = 36377 ; Dt = 0.0118s ; Time = 428.58s ; ETA = 64.21s Iter = 38109 ; Dt = 0.0118s ; Time = 448.98s ; ETA = 67.01s Iter = 39841 ; Dt = 0.0118s ; Time = 469.39s ; ETA = 56.65s Iter = 41574 ; Dt = 0.0118s ; Time = 489.81s ; ETA = 56.24s Iter = 43306 ; Dt = 0.0118s ; Time = 510.21s ; ETA = 53.07s Iter = 45038 ; Dt = 0.0118s ; Time = 530.62s ; ETA = 60.41s Iter = 46770 ; Dt = 0.0118s ; Time = 551.03s ; ETA = 48.25s Iter = 48502 ; Dt = 0.0118s ; Time = 571.43s ; ETA = 46.46s Iter = 50234 ; Dt = 0.0118s ; Time = 591.84s ; ETA = 44.11s Iter = 51966 ; Dt = 0.0118s ; Time = 612.25s ; ETA = 41.45s Iter = 53699 ; Dt = 0.0118s ; Time = 632.66s ; ETA = 39.80s Iter = 55431 ; Dt = 0.0118s ; Time = 653.07s ; ETA = 37.77s Iter = 57163 ; Dt = 0.0118s ; Time = 673.48s ; ETA = 41.30s Iter = 58895 ; Dt = 0.0118s ; Time = 693.88s ; ETA = 33.62s Iter = 60627 ; Dt = 0.0118s ; Time = 714.29s ; ETA = 31.62s Iter = 62359 ; Dt = 0.0118s ; Time = 734.70s ; ETA = 28.86s Iter = 64092 ; Dt = 0.0118s ; Time = 755.11s ; ETA = 26.28s Iter = 65824 ; Dt = 0.0118s ; Time = 775.52s ; ETA = 24.80s Iter = 67556 ; Dt = 0.0118s ; Time = 795.93s ; ETA = 23.04s Iter = 69288 ; Dt = 0.0118s ; Time = 816.33s ; ETA = 20.04s Iter = 71020 ; Dt = 0.0118s ; Time = 836.74s ; ETA = 17.60s Iter = 72752 ; Dt = 0.0118s ; Time = 857.15s ; ETA = 15.52s Iter = 74484 ; Dt = 0.0118s ; Time = 877.55s ; ETA = 13.15s Iter = 76217 ; Dt = 0.0118s ; Time = 897.97s ; ETA = 11.11s Iter = 77949 ; Dt = 0.0118s ; Time = 918.38s ; ETA = 8.97s Iter = 79681 ; Dt = 0.0118s ; Time = 938.78s ; ETA = 6.58s Iter = 81413 ; Dt = 0.0118s ; Time = 959.19s ; ETA = 4.38s Iter = 83145 ; Dt = 0.0118s ; Time = 979.60s ; ETA = 2.24s Iter = 84877 ; Dt = 0.0104s ; Time = 1000.00s ; ETA = 0.00s =================================================================== | END | =================================================================== Problem.solve() completed in 110.29126620292664s (wall time) SOLVING CASE 6 =================================================================== | INITIALIZATION | =================================================================== Problem size: Nodes=400, Layers=1, Triangles=634, Iter = 0 ; Dt = 0.0000s ; Time = 0.00s ; ETA = 0.00s =================================================================== | TIME LOOP | =================================================================== Iter = 1899 ; Dt = 0.0107s ; Time = 20.41s ; ETA = 31.18s Iter = 3797 ; Dt = 0.0108s ; Time = 40.82s ; ETA = 35.16s Iter = 5696 ; Dt = 0.0107s ; Time = 61.23s ; ETA = 35.45s Iter = 7594 ; Dt = 0.0108s ; Time = 81.64s ; ETA = 28.14s Iter = 9492 ; Dt = 0.0107s ; Time = 102.04s ; ETA = 30.43s Iter = 11391 ; Dt = 0.0108s ; Time = 122.46s ; ETA = 28.33s Iter = 13289 ; Dt = 0.0108s ; Time = 142.86s ; ETA = 25.91s Iter = 15187 ; Dt = 0.0108s ; Time = 163.27s ; ETA = 26.76s Iter = 17086 ; Dt = 0.0108s ; Time = 183.68s ; ETA = 24.71s Iter = 18984 ; Dt = 0.0108s ; Time = 204.09s ; ETA = 23.99s Iter = 20882 ; Dt = 0.0108s ; Time = 224.49s ; ETA = 24.20s Iter = 22780 ; Dt = 0.0108s ; Time = 244.90s ; ETA = 22.91s Iter = 24679 ; Dt = 0.0108s ; Time = 265.32s ; ETA = 23.33s Iter = 26577 ; Dt = 0.0108s ; Time = 285.72s ; ETA = 21.51s Iter = 28475 ; Dt = 0.0108s ; Time = 306.13s ; ETA = 21.96s Iter = 30373 ; Dt = 0.0108s ; Time = 326.53s ; ETA = 20.31s Iter = 32272 ; Dt = 0.0108s ; Time = 346.95s ; ETA = 19.90s Iter = 34170 ; Dt = 0.0108s ; Time = 367.35s ; ETA = 19.18s Iter = 36068 ; Dt = 0.0108s ; Time = 387.76s ; ETA = 20.69s Iter = 37967 ; Dt = 0.0108s ; Time = 408.17s ; ETA = 18.16s Iter = 39865 ; Dt = 0.0108s ; Time = 428.58s ; ETA = 17.26s Iter = 41763 ; Dt = 0.0108s ; Time = 448.98s ; ETA = 16.68s Iter = 43661 ; Dt = 0.0108s ; Time = 469.39s ; ETA = 17.76s Iter = 45560 ; Dt = 0.0108s ; Time = 489.81s ; ETA = 15.56s Iter = 47458 ; Dt = 0.0108s ; Time = 510.21s ; ETA = 15.60s Iter = 49356 ; Dt = 0.0108s ; Time = 530.62s ; ETA = 14.24s Iter = 51254 ; Dt = 0.0108s ; Time = 551.02s ; ETA = 14.00s Iter = 53153 ; Dt = 0.0108s ; Time = 571.44s ; ETA = 13.01s Iter = 55051 ; Dt = 0.0108s ; Time = 591.84s ; ETA = 12.33s Iter = 56949 ; Dt = 0.0108s ; Time = 612.25s ; ETA = 12.22s Iter = 58848 ; Dt = 0.0108s ; Time = 632.66s ; ETA = 11.19s Iter = 60746 ; Dt = 0.0108s ; Time = 653.07s ; ETA = 10.79s Iter = 62644 ; Dt = 0.0108s ; Time = 673.47s ; ETA = 9.94s Iter = 64542 ; Dt = 0.0108s ; Time = 693.88s ; ETA = 10.96s Iter = 66441 ; Dt = 0.0108s ; Time = 714.30s ; ETA = 8.67s Iter = 68339 ; Dt = 0.0108s ; Time = 734.70s ; ETA = 8.00s Iter = 70237 ; Dt = 0.0108s ; Time = 755.11s ; ETA = 7.45s Iter = 72135 ; Dt = 0.0108s ; Time = 775.51s ; ETA = 7.05s Iter = 74034 ; Dt = 0.0108s ; Time = 795.93s ; ETA = 6.15s Iter = 75932 ; Dt = 0.0108s ; Time = 816.33s ; ETA = 5.63s Iter = 77830 ; Dt = 0.0108s ; Time = 836.74s ; ETA = 5.10s Iter = 79729 ; Dt = 0.0108s ; Time = 857.15s ; ETA = 4.31s Iter = 81627 ; Dt = 0.0108s ; Time = 877.56s ; ETA = 3.71s Iter = 83525 ; Dt = 0.0108s ; Time = 897.96s ; ETA = 3.08s Iter = 85423 ; Dt = 0.0108s ; Time = 918.37s ; ETA = 2.58s Iter = 87322 ; Dt = 0.0108s ; Time = 938.78s ; ETA = 1.97s Iter = 89220 ; Dt = 0.0108s ; Time = 959.19s ; ETA = 1.28s Iter = 91118 ; Dt = 0.0108s ; Time = 979.60s ; ETA = 0.63s Iter = 93016 ; Dt = 0.0104s ; Time = 1000.00s ; ETA = 0.00s =================================================================== | END | =================================================================== Problem.solve() completed in 31.63927960395813s (wall time) SOLVING CASE 7 =================================================================== | INITIALIZATION | =================================================================== Problem size: Nodes=400, Layers=1, Triangles=634, Iter = 0 ; Dt = 0.0000s ; Time = 0.00s ; ETA = 0.00s =================================================================== | TIME LOOP | =================================================================== Iter = 1898 ; Dt = 0.0108s ; Time = 20.42s ; ETA = 135.99s Iter = 3794 ; Dt = 0.0108s ; Time = 40.83s ; ETA = 126.29s Iter = 5690 ; Dt = 0.0108s ; Time = 61.23s ; ETA = 125.11s Iter = 7586 ; Dt = 0.0108s ; Time = 81.64s ; ETA = 119.53s Iter = 9483 ; Dt = 0.0108s ; Time = 102.05s ; ETA = 121.29s Iter = 11379 ; Dt = 0.0108s ; Time = 122.46s ; ETA = 117.86s Iter = 13275 ; Dt = 0.0108s ; Time = 142.86s ; ETA = 110.75s Iter = 15171 ; Dt = 0.0108s ; Time = 163.27s ; ETA = 109.30s Iter = 17067 ; Dt = 0.0108s ; Time = 183.68s ; ETA = 107.40s Iter = 18963 ; Dt = 0.0108s ; Time = 204.09s ; ETA = 104.17s Iter = 20859 ; Dt = 0.0108s ; Time = 224.50s ; ETA = 102.06s Iter = 22755 ; Dt = 0.0108s ; Time = 244.90s ; ETA = 97.28s Iter = 24651 ; Dt = 0.0108s ; Time = 265.31s ; ETA = 95.91s Iter = 26547 ; Dt = 0.0108s ; Time = 285.72s ; ETA = 95.55s Iter = 28443 ; Dt = 0.0108s ; Time = 306.13s ; ETA = 90.48s Iter = 30339 ; Dt = 0.0108s ; Time = 326.54s ; ETA = 87.59s Iter = 32235 ; Dt = 0.0108s ; Time = 346.95s ; ETA = 85.11s Iter = 34131 ; Dt = 0.0108s ; Time = 367.36s ; ETA = 83.25s Iter = 36027 ; Dt = 0.0108s ; Time = 387.76s ; ETA = 80.32s Iter = 37923 ; Dt = 0.0108s ; Time = 408.17s ; ETA = 76.95s Iter = 39819 ; Dt = 0.0108s ; Time = 428.58s ; ETA = 75.45s Iter = 41715 ; Dt = 0.0108s ; Time = 448.99s ; ETA = 72.53s Iter = 43611 ; Dt = 0.0108s ; Time = 469.40s ; ETA = 68.73s Iter = 45507 ; Dt = 0.0108s ; Time = 489.81s ; ETA = 66.65s Iter = 47403 ; Dt = 0.0108s ; Time = 510.21s ; ETA = 67.47s Iter = 49298 ; Dt = 0.0108s ; Time = 530.61s ; ETA = 76.50s Iter = 51194 ; Dt = 0.0108s ; Time = 551.02s ; ETA = 59.28s Iter = 53090 ; Dt = 0.0108s ; Time = 571.43s ; ETA = 57.77s Iter = 54986 ; Dt = 0.0108s ; Time = 591.84s ; ETA = 58.57s Iter = 56882 ; Dt = 0.0108s ; Time = 612.25s ; ETA = 53.25s Iter = 58778 ; Dt = 0.0108s ; Time = 632.65s ; ETA = 48.29s Iter = 60674 ; Dt = 0.0108s ; Time = 653.06s ; ETA = 46.78s Iter = 62570 ; Dt = 0.0108s ; Time = 673.47s ; ETA = 42.97s Iter = 64466 ; Dt = 0.0108s ; Time = 693.88s ; ETA = 40.26s Iter = 66362 ; Dt = 0.0108s ; Time = 714.29s ; ETA = 37.84s Iter = 68258 ; Dt = 0.0108s ; Time = 734.70s ; ETA = 35.28s Iter = 70154 ; Dt = 0.0108s ; Time = 755.11s ; ETA = 32.44s Iter = 72050 ; Dt = 0.0108s ; Time = 775.51s ; ETA = 29.52s Iter = 73946 ; Dt = 0.0108s ; Time = 795.92s ; ETA = 26.85s Iter = 75842 ; Dt = 0.0108s ; Time = 816.33s ; ETA = 24.60s Iter = 77738 ; Dt = 0.0108s ; Time = 836.74s ; ETA = 21.82s Iter = 79634 ; Dt = 0.0108s ; Time = 857.15s ; ETA = 19.25s Iter = 81530 ; Dt = 0.0108s ; Time = 877.56s ; ETA = 16.28s Iter = 83426 ; Dt = 0.0108s ; Time = 897.96s ; ETA = 13.90s Iter = 85322 ; Dt = 0.0108s ; Time = 918.37s ; ETA = 10.90s Iter = 87218 ; Dt = 0.0108s ; Time = 938.78s ; ETA = 8.23s Iter = 89114 ; Dt = 0.0108s ; Time = 959.19s ; ETA = 5.35s Iter = 91010 ; Dt = 0.0108s ; Time = 979.60s ; ETA = 2.71s Iter = 92906 ; Dt = 0.0036s ; Time = 1000.00s ; ETA = 0.00s =================================================================== | END | =================================================================== Problem.solve() completed in 134.0907199382782s (wall time) .. GENERATED FROM PYTHON SOURCE LINES 220-224 Errors: -------------------- Errors computed during solving loop are stored in ``error_*`` matrix. .. GENERATED FROM PYTHON SOURCE LINES 225-235 .. code-block:: default error_H = np.zeros((N_mesh, 2)) error_QX = np.zeros((N_mesh, 2)) error_QY = np.zeros((N_mesh, 2)) for I, case in enumerate(cases): 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 .. GENERATED FROM PYTHON SOURCE LINES 236-237 Errors are save in h5 format in \data repertory. .. GENERATED FROM PYTHON SOURCE LINES 238-249 .. code-block:: default if SAVE_ERROR: os.system('rm data/bump_errors.h5') REF_FILE = "data/bump_errors.h5" with h5py.File(REF_FILE, 'w') as f: 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) .. GENERATED FROM PYTHON SOURCE LINES 250-253 Errors plots: -------------------- .. GENERATED FROM PYTHON SOURCE LINES 254-300 .. code-block:: default def plot_error(ax, cases, error): err = np.zeros((N_mesh)) ab = np.zeros((N_mesh)) for order in range(2): for M in range(N_mesh): err[M] = np.log(error[M, order]) ab[M] = 0.5*np.log(NT[M]/NT[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]), -1.0*xmax+np.log(error[0, 1])],\ 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"] = [14, 6] colors = ['blue','red','green','purple'] markers = ['s','*','o','v'] labels = ['$1^{st}$ order','$2^{nd}$ order'] fig = plt.figure() ax0 = fig.add_subplot(111) plot_error(ax0, cases, error_H) plot_reference(ax0, error_H) 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() fig = plt.figure() ax1 = fig.add_subplot(111) plot_error(ax1, cases, error_QX) plot_reference(ax1, error_QX) ax1.set_xlabel('$log(l_0/l_i)$') ax1.set_ylabel('Error of $hu$ in $L_2$ norm') ax1.set_title("Convergence of $hv$") plt.legend(loc=1) #plt.show() .. rst-class:: sphx-glr-horizontal * .. image-sg:: /auto_examples_mesh/images/sphx_glr_example_meshconvergence_bump_005.png :alt: Convergence of $h$ :srcset: /auto_examples_mesh/images/sphx_glr_example_meshconvergence_bump_005.png :class: sphx-glr-multi-img * .. image-sg:: /auto_examples_mesh/images/sphx_glr_example_meshconvergence_bump_006.png :alt: Convergence of $hv$ :srcset: /auto_examples_mesh/images/sphx_glr_example_meshconvergence_bump_006.png :class: sphx-glr-multi-img .. rst-class:: sphx-glr-timing **Total running time of the script:** ( 6 minutes 11.669 seconds) .. _sphx_glr_download_auto_examples_mesh_example_meshconvergence_bump.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_bump.py ` .. container:: sphx-glr-download sphx-glr-download-jupyter :download:`Download Jupyter notebook: example_meshconvergence_bump.ipynb ` .. only:: html .. rst-class:: sphx-glr-signature `Gallery generated by Sphinx-Gallery `_