Examples added, example readme and example unit test updated

This commit is contained in:
Nicolas Kruse 2025-06-06 13:15:34 +02:00
parent 0578b552be
commit 7671a89e14
8 changed files with 215 additions and 17 deletions

View File

@ -16,7 +16,7 @@ def run_cmd(command: list[str]):
assert (not result.stderr or assert (not result.stderr or
any('RuntimeWarning: ' in line for line in result.stderr.splitlines()) or any('RuntimeWarning: ' in line for line in result.stderr.splitlines()) or
any('[NbConvertApp]' in line for line in result.stderr.splitlines())), 'ERROR: ' + result.stderr any('[NbConvertApp]' in line and 'error' not in line.lower() for line in result.stderr.splitlines())), 'ERROR: ' + result.stderr
def run_rendering(input_path: str, output_directory: str): def run_rendering(input_path: str, output_directory: str):

View File

@ -21,5 +21,5 @@ jupyter nbconvert --to markdown docs/source/_autogenerated/soec_methane.ipynb --
A new example Markdown file can be created from a Jupyter Notebook running A new example Markdown file can be created from a Jupyter Notebook running
the following command: the following command:
``` bash ``` bash
jupyter nbconvert --to new_example.ipynb --NbConvertApp.use_output_suffix=False --ClearOutputPreprocessor.enabled=True --output-dir examples/ --output new_example.md jupyter nbconvert --to markdown new_example.ipynb --NbConvertApp.use_output_suffix=False --ClearOutputPreprocessor.enabled=True --output-dir examples/ --output new_example.md
``` ```

View File

@ -0,0 +1,68 @@
# Carbon Activity
This example shows the equilibrium calculation for solid carbon.
```python
import gaspype as gp
import numpy as np
import matplotlib.pyplot as plt
```
Setting temperatures and pressure:
```python
t_range = np.array([600, 700, 800, 900, 1100, 1500]) # °C
p = 1e5 # Pa
fs = gp.fluid_system(['H2', 'H2O', 'CO2', 'CO', 'CH4'])
```
Equilibrium calculation for methane steam mixtures:
```python
ratio = np.linspace(0.01, 1.5, num=128)
fl = gp.fluid({'CH4': 1}, fs) + ratio * gp.fluid({'H2O': 1}, fs)
```
gaspype.carbon_activity supports currently only 0D fluids therefore we build this helper function:
```python
def partial_c_activity(fl: gp.fluid, t: float, p: float):
fls = fl.array_composition.shape
eq_fl = gp.equilibrium(fl, t, p)
ret = np.zeros(fls[0])
for i in range(fls[0]):
ret[i] = gp.carbon_activity(gp.fluid(eq_fl.array_composition[i,:], fs), t, p)
return ret
```
Now we use the helper function to calculate the carbon activitx for all
compositions in equilibrium_h2o times all temperatures in t_range:
```python
carbon_activity = np.vstack([partial_c_activity(fl, tc + 273.15, p) for tc in t_range])
```
Plot carbon activities, a activity of > 1 means there is thermodynamically the formation of sold carbon favored.
```python
fig, ax = plt.subplots(figsize=(6, 4), dpi=120)
ax.set_xlabel("CO2/CH4")
ax.set_ylabel("carbon activity")
ax.set_ylim(1e-1, 1e3)
ax.set_yscale('log')
ax.plot(ratio, carbon_activity.T)
ax.hlines(1, np.min(ratio), np.max(ratio), colors='k', linestyles='dashed')
ax.legend([f'{tc} °C' for tc in t_range])
```

View File

@ -0,0 +1,60 @@
# Methane Mixtures
This example shows equilibria of methane mixed with steam and CO2
```python
import gaspype as gp
import numpy as np
import matplotlib.pyplot as plt
```
Setting temperature and pressure:
```python
t = 900 + 273.15
p = 1e5
fs = gp.fluid_system(['H2', 'H2O', 'CO2', 'CO', 'CH4', 'O2'])
```
Equilibrium calculation for methane steam mixtures:
```python
ratio = np.linspace(0.01, 1.5, num=64)
el = gp.fluid({'CH4': 1}, fs) + ratio * gp.fluid({'H2O': 1}, fs)
equilibrium_h2o = gp.equilibrium(el, t, p)
```
```python
fig, ax = plt.subplots(figsize=(6, 4), dpi=120)
ax.set_xlabel("H2O/CH4")
ax.set_ylabel("molar fraction")
ax.set_ylim(0, 1.1)
#ax.set_xlim(0, 100)
ax.plot(ratio, equilibrium_h2o.get_x())
ax.legend(fs.active_species)
```
Equilibrium calculation for methane CO2 mixtures:
```python
el = gp.fluid({'CH4': 1}, fs) + ratio * gp.fluid({'H2O': 1}, fs)
equilibrium_co2 = gp.equilibrium(el, t, p)
```
```python
fig, ax = plt.subplots(figsize=(6, 4), dpi=120)
ax.set_xlabel("CO2/CH4")
ax.set_ylabel("molar fraction")
ax.set_ylim(0, 1.1)
#ax.set_xlim(0, 100)
ax.plot(ratio, equilibrium_co2.get_x())
ax.legend(fs.active_species)
```

View File

@ -1,4 +1,8 @@
# SOEC example # SOEC with Methane
This example shows a 1D isothermal SOEC (Solid oxide electrolyzer cell) model.
The operating parameters chosen here are not necessary realistic
```python ```python
import gaspype as gp import gaspype as gp
@ -7,7 +11,8 @@ import numpy as np
import matplotlib.pyplot as plt import matplotlib.pyplot as plt
``` ```
Calculation of the local equilibrium compositions on the fuel and air
side in counter flow along the fuel flow direction:
```python ```python
fuel_utilization = 0.90 fuel_utilization = 0.90
air_utilization = 0.5 air_utilization = 0.5
@ -28,9 +33,8 @@ fuel_side = gp.equilibrium(feed_fuel + perm_oxygen, t, p)
air_side = gp.equilibrium(feed_air - perm_oxygen, t, p) air_side = gp.equilibrium(feed_air - perm_oxygen, t, p)
``` ```
Plot compositions of the fuel and air side:
```python ```python
#Plot compositions on fuel and air side
fig, ax = plt.subplots() fig, ax = plt.subplots()
ax.set_xlabel("Conversion") ax.set_xlabel("Conversion")
ax.set_ylabel("Molar fraction") ax.set_ylabel("Molar fraction")
@ -44,13 +48,14 @@ ax.plot(conversion, air_side.get_x(), '-')
ax.legend(air_side.species) ax.legend(air_side.species)
``` ```
Calculation of the oxygen partial pressures:
```python ```python
o2_fuel_side = gp.oxygen_partial_pressure(fuel_side, t, p) o2_fuel_side = gp.oxygen_partial_pressure(fuel_side, t, p)
o2_air_side = air_side.get_x('O2') * p o2_air_side = air_side.get_x('O2') * p
``` ```
Plot oxygen partial pressures:
```python ```python
#Plot oxygen partial pressure
fig, ax = plt.subplots() fig, ax = plt.subplots()
ax.set_xlabel("Conversion") ax.set_xlabel("Conversion")
ax.set_ylabel("Oxygen partial pressure / Pa") ax.set_ylabel("Oxygen partial pressure / Pa")
@ -59,13 +64,14 @@ ax.plot(conversion, np.stack([o2_fuel_side, o2_air_side], axis=1), '-')
ax.legend(['o2_fuel_side', 'o2_air_side']) ax.legend(['o2_fuel_side', 'o2_air_side'])
``` ```
Calculation of the local nernst potential between fuel and air side:
```python ```python
z_O2 = 4 z_O2 = 4
nernst_voltage = R*t / (z_O2*F) * np.log(o2_air_side/o2_fuel_side) nernst_voltage = R*t / (z_O2*F) * np.log(o2_air_side/o2_fuel_side)
``` ```
#Plot nernst potential:
```python ```python
#Plot voltage potential
fig, ax = plt.subplots() fig, ax = plt.subplots()
ax.set_xlabel("Conversion") ax.set_xlabel("Conversion")
ax.set_ylabel("Voltage / V") ax.set_ylabel("Voltage / V")
@ -73,6 +79,15 @@ ax.plot(conversion, nernst_voltage, '-')
print(np.min(nernst_voltage)) print(np.min(nernst_voltage))
``` ```
The model uses between each node a constant conversion. Because
current density depends strongly on the position along the cell
the constant conversion does not relate to a constant distance.
![Alt text](../../media/soc_inverted.svg)
To calculate the local current density (**node_current**) as well
as the total cell current (**terminal_current**) the (relative)
physical distance between the nodes (**dz**) must be calculated:
```python ```python
cell_voltage = 0.77 #V cell_voltage = 0.77 #V
ASR = 0.2 #Ohm*cm² ASR = 0.2 #Ohm*cm²
@ -88,8 +103,8 @@ terminal_current = np.sum(current * dz) #mA/cm² (Total cell current per cell ar
print(f'Terminal current: {terminal_current:.2f} A/cm²') print(f'Terminal current: {terminal_current:.2f} A/cm²')
``` ```
Plot the local current density:
```python ```python
#Plot current density
z_position = np.concatenate([[0], np.cumsum(dz)]) # Relative position of each node z_position = np.concatenate([[0], np.cumsum(dz)]) # Relative position of each node
fig, ax = plt.subplots() fig, ax = plt.subplots()

View File

@ -0,0 +1,52 @@
# Sulfur Oxygen Equilibrium
This example shows equilibrium calculations for sulfur/oxygen mixtures.
```python
import gaspype as gp
import numpy as np
import matplotlib.pyplot as plt
```
List possible sulfur/oxygen species:
```python
gp.species(element_names = 'S, O')
```
Or more specific by using regular expressions:
```python
gp.species('S?[2-3]?O?[2-5]?', use_regex=True)
```
Calculation of the molar equilibrium fractions for sulfur and oxygen depending on the oxygen to sulfur ratio:
```python
fs = gp.fluid_system(['S2', 'S2O', 'SO2', 'SO3', 'O2'])
oxygen_ratio = np.linspace(0.5, 3, num=128)
el = gp.elements({'S': 1}, fs) + oxygen_ratio * gp.elements({'O': 1}, fs)
composition = gp.equilibrium(el, 800+273.15, 1e4)
plt.plot(oxygen_ratio, composition.get_x())
plt.legend(composition.species)
```
Calculation of the molar equilibrium fractions for sulfur and oxygen depending on temperature in °C:
```python
fs = gp.fluid_system(['S2', 'S2O', 'SO2', 'SO3', 'O2'])
el = gp.elements({'S': 1, 'O':2.5}, fs)
t_range = np.linspace(500, 1300, num=32)
composition = gp.equilibrium(el, t_range+273.15, 1e4)
plt.plot(t_range, composition.get_x())
plt.legend(composition.species)
```

View File

@ -93,7 +93,10 @@ def segments_to_test(segments: Iterable[markdown_segment], script_language: str
if segment.code_block: if segment.code_block:
if segment.language == script_language: if segment.language == script_language:
lines = [line for line in segment.text.splitlines() if line.strip()] lines = [line for line in segment.text.splitlines() if line.strip()]
ret_block_flag = lines[-1] if not re.match(r'^[^(]*=', lines[-1]) and not lines[-1].startswith('import ') else None ret_block_flag = lines[-1] if (not re.match(r'^[^(]*=', lines[-1]) and
not lines[-1].startswith('import ') and
not lines[-1].startswith('from ') and
not lines[-1].startswith(' ')) else None
# print('Last line: ', ret_block_flag, '-----------', lines[-1]) # print('Last line: ', ret_block_flag, '-----------', lines[-1])
yield '' yield ''

View File

@ -11,7 +11,7 @@ def test_readme():
def test_example_code(): def test_example_code():
filter = 'docs/source/examples/*.md' filter = 'examples/*.md'
files = glob(filter) files = glob(filter)
for path in files: for path in files:
@ -19,7 +19,7 @@ def test_example_code():
if not file_name.lower() == 'readme': if not file_name.lower() == 'readme':
print(f"> Test Example {file_name} ...") print(f"> Test Example {file_name} ...")
md_to_code.convert_to('test', path, f'tests/autogenerated_{file_name}.py') md_to_code.convert_to('test', path, f'tests/autogenerated_{file_name}.py')
mod = importlib.import_module(file_name) mod = importlib.import_module(f'autogenerated_{file_name}')
mod.run_test() mod.run_test()