Workflow: IBOFlow

[anderslogg]

Workflow for supporting simulations with IPS IBOFlow. See workflows/workflow_iboflow.py in DTCC module.

Questions to be answered:

Overview: Explain what the overall idea is

Input data: List input data, file formats etc in detail. How to obtain the data? Also attach example data to this issue.

Output data: List input data, file formats etc in detail. How is the output data used?

How is this done today: Explain the workflow as it is today with existing tools.

Comments: Any comments? How would this work in an ideal world?

[anderslogg]

Input from Franziska:

Input data

XML Setup file, VTK or CityJSON files for terrain and building triangulations, windrose statistics from nearby weather station (optionally), additional building or terrain triangulations for exploring other design options (optionally)

Output data

XDMF + HDF5 files for the entire discretized flow field, surface meshes, streamlines, point data

How is this done today

XML file generated manually, default settings are set so that a minimum of input settings are required by the user (We can provide a minimum exemplary XML input file).

VTK or CityJSON is usually generated via DTCC Builder, geometries cropped internally.

Windrose is usually generated by processing the data downloaded for a chosen weather station from SMHI (I think, coupling the DTCC Platform to SMHI for direct download or data visualization would be super useful as very many data can be provided that is of interest within the DTC) using a python script which can be provided to the DTCC.

Simulations are either run on the local desktop computer or send to the queue on the computing cluster.

Vis in paraview and UE.

How would this work in an ideal world

In a GUI mark the area that you are interested in and therewith automatically download the geometry data from Lantmäteriet and process it for generating surface triangulation.

Automatically grab the marked area as simulation domain area

Let the user specify in the GUI the minimum input, like discretization, type of simulation (i.e. single run or wind comfort).

Let the user specify which weather station location to choose from and which time interval to take and provide the windrose statistics.

Run the simulation.

Vis the simulation results.

[anderslogg]

Input on getting weather data:

Choose measurement station from https://www.smhi.se/data/meteorologi/ladda-ner-meteorologiska-observationer/#param=wind,stations=core

Download data contained in: Ladda ner historisk data, förutom de senaste 3 månaderna (.csv)

Process the data with python script to obtain a windrose file

#!/usr/bin/env python3
​
import numpy as np
import pandas as pd
import math
import datetime
import matplotlib.pyplot as plt
import argparse
​
from windrose import WindroseAxes
​
​
""" Script to generate the windrose from measruements of wind speeds according to the 
SMHI format:
Datum       Tid (UTC)   Vindriktning    Kvalitet    Vindhastighet   Kvalitet
1951-01-01  00:00:00    100                 G           3.0             G
1951-01-01  03:00:00    110                 G           2.0             G
1951-01-01  06:00:00    70                  G           2.0             G
​
"""
​
parser = argparse.ArgumentParser(
                    prog='CreateWindrose',
                    description='Creates windrose statistics from SMHI data file in xml format for IBOFlow and plot in addition',
                    epilog='No warranty')
​
parser.add_argument('filename')           # positional argument
parser.add_argument('-nDeg', '--numberDegrees')      # option that takes a value
​
args = parser.parse_args()
​
class opts:
    """ Class containing the options"""
    excludeDefect=False #needs always to be set to false for now
    nDeg = int(args.numberDegrees)
    nVel = 6
​
#Read and convert the input data
#with open("smhi-opendata_3_4_71420_20220111_074056.csv") as file:
#    data = np.loadtxt(file, skiprows=25, usecols=(0,1,2,3,4,5), delimiter = ";", dtype=str, max_rows=23000)
data = pd.read_csv(args.filename, skiprows=15, delimiter=";", usecols=(0,1,2,3,4,5))
print(data)
​
height = 10.0
date = data['Datum']
time = data['Tid (UTC)']
timeObj = np.zeros(len(time),dtype=datetime.datetime)
k=0
for d,t in zip(date,time):
    timeStr = d+" "+t
    timeObj[k] = datetime.datetime.strptime(timeStr, '%Y-%m-%d %H:%M:%S')
    k=k+1
direction = data['Vindriktning'].astype(float) 
qualityD = data['Kvalitet'] #G for good and y for defective
speed = data['Vindhastighet'].astype(float)
qualityS = data['Kvalitet.1'] # G for good and y for defective
​
#Create the binning
binSize=360/opts.nDeg
binDeg=np.arange(0,360+binSize,binSize)
​
#binVelSize=(max(speed)-min(speed))/(opts.nVel-1)
binVel=[0, 0.5, 2, 4, 6, 8, 10, 12, 14, 16, 18, 100]
opts.nVel=len(binVel)-1
#binVel=np.arange(min(speed)-max(binVelSize,0),max(speed)+binVelSize,binVelSize)
windrose = np.zeros((opts.nDeg,opts.nVel))
​
print(binVel)
print(binDeg)
​
# Put all the data in the bins
k=0
totalTime=0
count=0
for d,qD,s,qS,t in zip(direction,qualityD,speed,qualityS,timeObj):
    if not (opts.excludeDefect and (qD=='Y' or qS=='Y')): #TODO: if excludeDefect is used, the time needs to adjusted
        #print(s)
        i=0
        success=False
        while i<opts.nDeg:
            if math.isclose(d,0):
                d=1 #Hack to assign 0 also to the smallest bin
            if (d<=binDeg[i+1] and d>binDeg[i]):
                j=0
                while j<opts.nVel:
                    if ((s<=binVel[j+1] and s>binVel[j]) or (math.isclose(s,binVel[j+1]))):
                        success=True
                        #print(s)
                        if t.hour==0 or k==0:
                            windrose[i,j]=windrose[i,j]+(timeObj[k+1].hour-t.hour) #Assuming full hours
                            totalTime=totalTime+(timeObj[k+1].hour-t.hour)
                        elif k==len(direction) and t.hour==0:
                            pass #Do nothing
                        else:
                            windrose[i,j]=windrose[i,j]+(t.hour-timeObj[k-1].hour) #Assuming full hours
                            totalTime=totalTime+(t.hour-timeObj[k-1].hour)
                    j=j+1
            i=i+1
    if success: count=count+1
    k=k+1
print("Successfull findings: ",count)
​
#Calculate the frequency
print("Total time:", totalTime)
print("Windrose:", windrose)
print("Total time:", sum(sum(windrose)))
print("k:", k)
windrose=windrose/totalTime
#
​
f = open("windrose.xml", "w")
f.write("<windrose_data>\n\t<height value =\" ")
f.write(str(height))
f.write("\" /> \n")
f.write("\t<reference_roughness value = \"0.5\" \> \n")
f.write("\t<roughness_list>\n")
l=1
while l<len(binDeg): #binning for number of bins -1 or choose the binning mean
    f.write("\t\t<roughness value = \"")
    f.write(str(0.5))
    f.write("\"/> \n ")
    l=l+1
f.write("\t</roughness_list>\n")
f.write("\t<angle_list>\n")
l=1
while l<len(binDeg): #binning for number of bins -1 or choose the binning mean
    f.write("\t\t<angle value = \"")
    f.write(str(binDeg[l]))
    f.write("\"/> \n ")
    l=l+1
f.write("\t</angle_list>\n")
f.write("\t<velocity_list>\n")
l=1
while l<len(binVel): #binning for number of bins -1 or choose the binning mean
    f.write("\t\t<velocity x = \"")
    f.write(str(round(binVel[l],2)))
    f.write("\" y = \"0.0\" z = \"0.0\" /> \n ")
    l=l+1
f.write("\t</velocity_list>\n")
k=0
f.write("\t<frequency>\n")
while k<windrose.shape[0]:
    f.write("\t\t<angle>\n")
    l=0
    while l<windrose.shape[1]:
        f.write("\t\t\t<vel value = \"")
        f.write(str(round(windrose[k][l],5)))
        f.write("\"/> \n ")
        l=l+1
    f.write("\t\t</angle>\n")
    k=k+1
f.write("\t</frequency>\n")
f.write("</windrose_data>\n")
f.close()
​
f=plt.figure
ax=WindroseAxes.from_ax()
ax.bar(direction,speed,nsector=12)
ax.legend()
plt.savefig("Windrose.png",format='png')