KVSlab / VaMPy

A collection of tools for pre-processing, simulating, and post-processing vascular morphologies.
https://kvslab.github.io/VaMPy
GNU General Public License v3.0
15 stars 8 forks source link

Issue with pre-boundary_conditions with Oasis+VaMPy #141

Closed perezrmaria closed 6 months ago

perezrmaria commented 8 months ago

Hi! We're a team trying to run a simulation on the left atrium and have encountered an error that we're not able to solve.

We've used VaMPy to preprocess our mesh, and previously tried with a mesh provided as a VaMPy test (as we assume that mesh is okay) and after running oasis NSfracStep problem=Artery mesh_path=models/artery/artery.xml.gz save_solution_after_cycle=0 we've run into this error:

=== Mesh information === X range: 25.1176 to 42.5953 (delta: 17.4777) Y range: 24.3107 to 37.2822 (delta: 12.9715) Z range: 22.3903 to 44.4398 (delta: 22.0495) Number of cells: 158380 Number of cells per processor: 158380 Number of edges: 0 Number of faces: 320834 Number of facets: 320834 Number of vertices: 27982 Volume: 249.2360 Number of cells per volume: 635.4619 Creating initial folders === Initial pressure and area fraction === Boundary ID=2, pressure: 0.32641, area fraction: 0.33467 Boundary ID=3, pressure: 0.67359, area fraction: 0.66533 Traceback (most recent call last): File "/Users/maria/miniconda3_86/envs/vampy_env/bin/oasis", line 8, in <module> sys.exit(main()) ^^^^^^ File "/Users/maria/miniconda3_86/envs/vampy_env/lib/python3.11/site-packages/oasis/run_oasis.py", line 11, in main from oasis import NSfracStep File "/Users/maria/miniconda3_86/envs/vampy_env/lib/python3.11/site-packages/oasis/NSfracStep.py", line 172, in <module> vars().update(pre_solve_hook(**vars())) ^^^^^^^^^^^^^^^^^^^^^^^^ File "/Users/maria/miniconda3_86/envs/vampy_env/lib/python3.11/site-packages/oasis/problems/NSfracStep/Artery.py", line 177, in pre_solve_hook eval_dict["centerline_u_x_probes"] = Probes(probe_points.flatten(), V) ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^ TypeError: Probes() takes no arguments

We've tried with our own mesh and got the exact same error, for which we haven't found anything in the documentation. If someone has experienced this error before or has any knowledge on how to fix it we would be really happy to use some advice!

Thank you so much!

keiyamamo commented 8 months ago

hi @perezrmaria

There must be something wrong when importing Probes class at the beginning of the problem file because Probes should take arguments if they are imported properly. Have you seen something like this in the terminal Failed to import probe.probe11 ?

perezrmaria commented 8 months ago

Hi!

Yes, we've encountered that error and it says exactly that `Importing problem module Artery: /Users/maria/miniconda3_86/envs/vampy_env/lib/python3.11/site-packages/oasis/problems/NSfracStep/Artery.py

Failed to import probe.probe11

Start 274 MB 274 MB 36449 MB 36449 MB === Starting simulation for case: artery ===`

We're not sure about how we can import that module, as if we call conda list (since we've installed VaMPy through conda, that module is not installed. We have the probe folder, but are not sure on how to use it or import it. We've tried to look for it on the internet and the only documentation we found is this link (https://pypi.org/project/probe/) but it's outdated and we've not been able to install it (we don't think this is what we should use).

Any tips on how to solve this issue? Thank you so much!

keiyamamo commented 8 months ago

Probe class is part of vampy and is located here https://github.com/KVSlab/VaMPy/tree/master/src/vampy/simulation/probe they are written in c++ and not python, so it can cause problem sometime. Unfortunately I don’t know exactly how to fix the issue immediately because I didn’t encounter the problem myself. What you can do is to simply comment out all the lines related to probes. Probes are simply evaluating the velocity and pressure at points during the simulation and does not affect the CFD solution. It is always possible to get the point values after simulation has completed.

perezrmaria commented 8 months ago

Hi,

Thanks for the idea, we've commented all lines related to probes and obtained results for the artery problem. But from this, some doubts have arisen:

Any advice will be appreciated, thank you so much!

keiyamamo commented 8 months ago

Lines related to probes should not affect the CFD since all the CFD computations are done inside the solver, which is more or less independent from the problem file itself. I have two questions.

  1. what was the flow rate printed during the simulation? when running artery problem, you should see flow rate at the inlet. Do they make sense? If those values are not logical, something must be wrong in your problem file.
  2. can you attach problem file and visualization files for the velocity ?
keiyamamo commented 8 months ago

For your own mesh, can you provide more info about the mesh? can you copy paste the mesh information shown by vampy?

perezrmaria commented 8 months ago

Hi! Thank you so much for your response!

  1. During the simulation we've obtained flow rates that didn't look abnormal, but we don't have a lot of experience in this field. Here's a little peak: ` ========== Time step 230 ========== Sum of Q_out = 1.1764, Q_in = 1.1820, mean velocity (inlet): 0.1599, Reynolds number (inlet): 148.5580 For outlet with boundary ID=2: target flow rate: 0.3956 mL/s, computed flow rate: 0.4037 mL/s, pressure updated to: 0.2592 For outlet with boundary ID=3: target flow rate: 0.7864 mL/s, computed flow rate: 0.7727 mL/s, pressure updated to: 0.2855

========== Time step 240 ========== Sum of Q_out = 1.1998, Q_in = 1.2055, mean velocity (inlet): 0.1630, Reynolds number (inlet): 151.5093 For outlet with boundary ID=2: target flow rate: 0.4034 mL/s, computed flow rate: 0.4023 mL/s, pressure updated to: 0.2859 For outlet with boundary ID=3: target flow rate: 0.8021 mL/s, computed flow rate: 0.7975 mL/s, pressure updated to: 0.2537

========== Time step 250 ========== Sum of Q_out = 1.2239, Q_in = 1.2297, mean velocity (inlet): 0.1663, Reynolds number (inlet): 154.5468 For outlet with boundary ID=2: target flow rate: 0.4115 mL/s, computed flow rate: 0.4009 mL/s, pressure updated to: 0.2395 For outlet with boundary ID=3: target flow rate: 0.8181 mL/s, computed flow rate: 0.8230 mL/s, pressure updated to: 0.2580

========== Time step 260 ========== Sum of Q_out = 1.2487, Q_in = 1.2545, mean velocity (inlet): 0.1697, Reynolds number (inlet): 157.6696 For outlet with boundary ID=2: target flow rate: 0.4199 mL/s, computed flow rate: 0.4114 mL/s, pressure updated to: 0.1854 For outlet with boundary ID=3: target flow rate: 0.8347 mL/s, computed flow rate: 0.8373 mL/s, pressure updated to: 0.2730

========== Time step 270 ========== Sum of Q_out = 1.2741, Q_in = 1.2800, mean velocity (inlet): 0.1731, Reynolds number (inlet): 160.8767 For outlet with boundary ID=2: target flow rate: 0.4284 mL/s, computed flow rate: 0.4285 mL/s, pressure updated to: 0.1690 For outlet with boundary ID=3: target flow rate: 0.8516 mL/s, computed flow rate: 0.8456 mL/s, pressure updated to: 0.2665

========== Time step 280 ========== Sum of Q_out = 1.3002, Q_in = 1.3062, mean velocity (inlet): 0.1767, Reynolds number (inlet): 164.1670 For outlet with boundary ID=2: target flow rate: 0.4372 mL/s, computed flow rate: 0.4421 mL/s, pressure updated to: 0.1829 For outlet with boundary ID=3: target flow rate: 0.8691 mL/s, computed flow rate: 0.8580 mL/s, pressure updated to: 0.2392

========== Time step 290 ========== Sum of Q_out = 1.3269, Q_in = 1.3331, mean velocity (inlet): 0.1803, Reynolds number (inlet): 167.5392 For outlet with boundary ID=2: target flow rate: 0.4461 mL/s, computed flow rate: 0.4481 mL/s, pressure updated to: 0.2000 For outlet with boundary ID=3: target flow rate: 0.8869 mL/s, computed flow rate: 0.8788 mL/s, pressure updated to: 0.2136

========== Time step 300 ========== Sum of Q_out = 1.3543, Q_in = 1.3605, mean velocity (inlet): 0.1840, Reynolds number (inlet): 170.9921 For outlet with boundary ID=2: target flow rate: 0.4553 mL/s, computed flow rate: 0.4503 mL/s, pressure updated to: 0.1919 For outlet with boundary ID=3: target flow rate: 0.9052 mL/s, computed flow rate: 0.9040 mL/s, pressure updated to: 0.2035

Time = 2.8530e+01, timestep = 300, End time = 1.9020e+03 Total computing time on previous 100 timesteps = 92.233040 `

  1. We're attaching the problem file in a zip file for the artery. We've been running it inside of oasis/problems/NsFracStep Artery.py.zip

  2. On the information about our mesh, this is the information the program shows when trying to run the simulation: === Mesh information === X range: -101.262 to 61.0634 (delta: 162.3254) Y range: 56.0725 to 144.566 (delta: 88.4935) Z range: -1313.45 to -1148.44 (delta: 165.0100) Number of cells: 2006583 Number of cells per processor: 2006583 Number of edges: 0 Number of faces: 4040701 Number of facets: 4040701 Number of vertices: 335798 Volume: 220649.6306 Number of cells per volume: 9.0940

keiyamamo commented 8 months ago

Based on the log file, simulation seems to go well. Could you attach xdmf and h5 files for the velocity? Remember you need to select right quantity in Paraview to visualize the result.

Thanks for sharing the information about the mesh. This looks good, but it is hard to tell what’s the problem without further information. Maybe could you attach the mesh here?

perezrmaria commented 8 months ago

Sure, thank you for you help. We attach the xmdf and h5 files for both velocity and pressure as well as the mesh and json files obtained from vampy-mesh in this link to Google Drive as they were too big to attach them here. https://drive.google.com/drive/folders/1LipNk7QwQp1opQsUZIddgu5hAurW7BQx?usp=sharing

We've realized that the info for inlets and outlets info was upside down, as it considered 4 outlets and 1 inlet, instead of seeing the 4 PVs as inlets and the MV as outlet. We manually tried to correct that in the json.

keiyamamo commented 8 months ago

The velocity looks correct. Please see the attached image. Make sure to change the coloring scheme from solid color to the quantity you want to visualize. (red box)

Screenshot 2024-03-03 at 15 15 53

When meshing the atrium, you need to use the flag -at in the command line to automatically switch inlet and outlet. https://github.com/KVSlab/VaMPy/blob/6c6468a4d8fcda89c330c9b01cea253e1bd3c8c3/src/vampy/automatedPreprocessing/automated_preprocessing.py#L525

perezrmaria commented 8 months ago

ups! Rockie error, sorry we should have seen that on the velocity! Thanks for the tip on the inlet/oulet switch. We've continued trying to simulate and this error keeps appearing, we don't know if we should try with a simpler mesh with less elements or if it's a problem of the code: Start simulations Traceback (most recent call last): File "/Users/maria/miniconda3_86/envs/vampy_env/bin/oasis", line 8, in <module> sys.exit(main()) ^^^^^^ File "/Users/maria/miniconda3_86/envs/vampy_env/lib/python3.11/site-packages/oasis/run_oasis.py", line 11, in main from oasis import NSfracStep File "/Users/maria/miniconda3_86/envs/vampy_env/lib/python3.11/site-packages/oasis/NSfracStep.py", line 200, in <module> velocity_tentative_solve(**vars()) File "/Users/maria/miniconda3_86/envs/vampy_env/lib/python3.11/site-packages/oasis/solvers/NSfracStep/IPCS_ABCN.py", line 230, in velocity_tentative_solve u_sol.solve(A, x_[ui], b[ui]) RuntimeError: Thank you so much for your time and patience!!

keiyamamo commented 8 months ago

Yes, I highly recommend starting from small number of elements and simple geometry if possible. Using simple geometry makes every small step faster and will be easy for us to. The error can come from many sources. Right now, the error can come from many sources and is hard to pin point.

keiyamamo commented 7 months ago

Should be fixed with #142 ?

hkjeldsberg commented 6 months ago

Should be fixed with #142 !