I am new to the field of FCD and paraview, so bare with me
I ran a wind simulation for an urban size for single wind direction (using intermediary Eddy3D in Rhino3D for BlueCFD).
The question is how can I visualize the pedestrian wind‐comfort (Nen8100) instead of raw m/s. I couldn’t find any resources by other how this can be done, so i try to run this code in python shell, but that doesn’t seem to work. Any advice is more then welcome. Nice to join this awesome community. I am on v 5.6.2 on windows
#!/usr/bin/env pvpython
'''
comfort_classification_nen8100.py
Usage:
pvpython comfort_classification_nen8100.py /path/to/case.foam
or in ParaView GUI: Tools -> Python Shell -> exec(open('comfort_classification_nen8100.py').read())
This script:
1. Loads a .foam case file with multiple time steps
2. Computes velocity magnitude (Umag)
3. Generates a binary indicator for exceedance of 5 m/s at each timestep
4. Uses TemporalStatistics to compute the fraction of timesteps exceeding 5 m/s
5. Classifies points into NEN 8100 comfort classes based on exceedance probability
6. Sets up an annotated color map and renders the result
'''
import sys
from paraview.simple import *
# 1. Load the OpenFOAM case with all time steps
foam_file = sys.argv[1] if len(sys.argv) > 1 else 'case.foam'
foam = OpenFOAMReader(FileName=[foam_file])
foam.CellArrays = ['U']
foam.SkipZeroTime = False
foam.UpdatePipeline()
# 2. Compute velocity magnitude
calc1 = Calculator(Input=foam)
calc1.ResultArrayName = 'Umag'
calc1.Function = 'mag(U)'
calc1.UpdatePipeline()
# 3. Compute exceedance indicator (1 if Umag > 5 m/s, else 0)
calc2 = Calculator(Input=calc1)
calc2.ResultArrayName = 'exceed5'
calc2.Function = 'iif(Umag > 5.0, 1.0, 0.0)'
calc2.UpdatePipeline()
# 4. Compute temporal statistics (mean of exceedance gives probability)
ts = TemporalStatistics(Input=calc2)
ts.UpdatePipeline()
# 5. Classify based on exceedance probability (Mean_exceed5 * 100)
pf = ProgrammableFilter(Input=ts)
pf.Script = '''
import numpy as np
from paraview.util.numpy_support import vtk_to_numpy, numpy_to_vtk
# Mean_exceed5 is fraction of time steps where Umag > 5
p = vtk_to_numpy(inputs[0].PointData['Mean_exceed5']) * 100.0
cat = np.zeros_like(p, dtype=np.int8)
# NEN 8100 classes:
# 0: A (< 2.5% exceedance)
# 1: B (2.5–5% exceedance)
# 2: C (5–10% exceedance)
# 3: D (10–20% exceedance)
# 4: E (>= 20% exceedance)
cat[p < 2.5] = 0
cat[(p >= 2.5) & (p < 5.0)] = 1
cat[(p >= 5.0) & (p < 10.0)] = 2
cat[(p >= 10.0) & (p < 20.0)] = 3
cat[p >= 20.0] = 4
vtk_arr = numpy_to_vtk(cat)
vtk_arr.SetName('ComfortClass')
output.PointData.AddArray(vtk_arr)
'''
pf.UpdatePipeline()
# 6. Visualize the classified result
view = GetActiveViewOrCreate('RenderView')
disp = Show(pf, view)
ColorBy(disp, ('POINTS', 'ComfortClass'))
# 7. Setup an annotated indexed color map
lut = GetColorTransferFunction('ComfortClass')
lut.IndexedLookup = 1
lut.Annotations = [
'0', 'A: < 2.5% exceedance',
'1', 'B: 2.5–5% exceedance',
'2', 'C: 5–10% exceedance',
'3', 'D: 10–20% exceedance',
'4', 'E: ≥ 20% exceedance'
]
lut.RescaleTransferFunction(0, 4)
lut.ShowAnnotations = 1
view.ResetCamera()
Render()