From 53e731a3fdaafe8305a3d3a8038e12e41b4054da Mon Sep 17 00:00:00 2001 From: Garrett Barter Date: Thu, 1 Jul 2021 12:49:32 -0600 Subject: [PATCH] Squashed 'RAFT/' changes from 96936f0f3..c32c2adca c32c2adca Output OpenFAST-compatible .1 and .3 files from HAMS: 31ce20d90 Merge pull request #5 from WISDEM/omdao_gb git-subtree-dir: RAFT git-subtree-split: c32c2adca607e2f25ce404cf3c2525bf87e69322 --- designs/OC3spar.yaml | 2 + raft/raft_fowt.py | 152 +++++++++++++++++++++++++++++++------------ raft/raft_model.py | 27 +++++++- 3 files changed, 137 insertions(+), 44 deletions(-) diff --git a/designs/OC3spar.yaml b/designs/OC3spar.yaml index 81a661830..2bdcc59ba 100644 --- a/designs/OC3spar.yaml +++ b/designs/OC3spar.yaml @@ -1064,6 +1064,8 @@ turbine: platform: + min_freq_BEM : 0.03 # [Hz] lowest frequency and frequency interval to use in BEM analysis + members: # list all members here - name : center_spar # [-] an identifier (no longer has to be number) diff --git a/raft/raft_fowt.py b/raft/raft_fowt.py index c3a63cfcf..483932e3e 100644 --- a/raft/raft_fowt.py +++ b/raft/raft_fowt.py @@ -2,6 +2,7 @@ import os import numpy as np +from scipy.interpolate import interp1d import pyhams.pyhams as ph import raft.member2pnl as pnl @@ -17,7 +18,7 @@ class FOWT(): '''This class comprises the frequency domain model of a single floating wind turbine''' - def __init__(self, design, w=[], mpb=None, depth=600): + def __init__(self, design, w, mpb, depth=600): '''This initializes the FOWT object which contains everything for a single turbine's frequency-domain dynamics. The initializiation sets up the design description. @@ -26,22 +27,24 @@ def __init__(self, design, w=[], mpb=None, depth=600): design : dict Dictionary of the design... w - Array of frequencies to be used in analysis + Array of frequencies to be used in analysis (rad/s) + mpb + A MoorPy Body object that represents this FOWT in MoorPy + depth + Water depth, positive-down. (m) ''' # basic setup self.nDOF = 6 self.Xi0 = np.zeros(6) # mean offsets of platform, initialized at zero [m, rad] - if len(w)==0: - w = np.arange(.01, 3, 0.01) # angular frequencies tp analyze (rad/s) - self.depth = depth self.w = np.array(w) - self.nw = len(w) # number of frequencies + self.nw = len(w) # number of frequencies + self.dw = w[1]-w[0] # frequency increment [rad/s] - self.k = np.array( [ waveNumber(w, self.depth) for w in self.w] ) # wave number + self.k = np.array([waveNumber(w, self.depth) for w in self.w]) # wave number [m/rad] self.rho_water = getFromDict(design['site'], 'rho_water', default=1025.0) @@ -50,13 +53,17 @@ def __init__(self, design, w=[], mpb=None, depth=600): design['turbine']['tower']['dlsMax'] = getFromDict(design['turbine']['tower'], 'dlsMax', default=5.0) + + potModMaster = getFromDict(design['platform'], 'potModMaster', dtype=int, default=0) + dlsMax = getFromDict(design['platform'], 'dlsMax' , default=5.0) + min_freq_BEM = getFromDict(design['platform'], 'min_freq_BEM', default=self.dw/2/np.pi) + self.dw_BEM = 2.0*np.pi*min_freq_BEM + #self.pyHAMS_use_radps = getFromDict(design['platform'], 'pyHAMS_use_radps', dtype=bool, default=False) # could support this option if needed + # member-based platform description self.memberList = [] # list of member objects - potModMaster = getFromDict(design['platform'], 'potModMaster', dtype=int, default=0) - dlsMax = getFromDict(design['platform'], 'dlsMax', default=5.0) - for mi in design['platform']['members']: if potModMaster==1: @@ -341,9 +348,27 @@ def calcStatics(self): - def calcBEM(self): - '''This generates a mesh for the platform and runs a BEM analysis on it. - The mesh is only for non-interesecting members flagged with potMod=1.''' + def calcBEM(self, FAST_outname='', dw=0, wMax=0, wInf=10.0): + '''This generates a mesh for the platform and runs a BEM analysis on it + using pyHAMS. It can also write adjusted .1 and .3 output files suitable + for use with OpenFAST. + The mesh is only made for non-interesecting members flagged with potMod=1. + + PARAMETERS + ---------- + FAST_outname : string + If not empty, OpenFAST/WAMIT-style .1 and .3 files will be written. + The provided string should be the full path and file name (minus extension) for + these .1 and .3 files to be written to. + dw : float + Optional specification of custom frequency increment (rad/s). + wMax : float + Optional specification of maximum frequency for BEM analysis (rad/s). Will only be + used if it is greater than the maximum frequency used in RAFT. + wInf : float + Optional specification of large frequency to use as approximation for infinite + frequency in pyHAMS analysis (rad/s). + ''' # desired panel size (longitudinal and azimuthal) dz = 3 @@ -375,45 +400,88 @@ def calcBEM(self): pnl.writeMeshToGDF(vertices) # also a GDF for visualization - # TODO: maybe create a 'HAMS Project' class: - # - methods: - # - create HAMS project structure - # - write HAMS input files - # - call HAMS.exe - # - attributes: - # - addedMass, damping, fEx coefficients - ph.create_hams_dirs(meshDir) + ph.create_hams_dirs(meshDir) # + + ph.write_hydrostatic_file(meshDir) # HAMS needs a hydrostatics file, but it's unused for .1 and .3, so write a blank one + + # prepare frequency settings for HAMS + if dw == 0: + dw_HAMS = self.dw_BEM + else: + dw_HAMS = dw # allow override of frequency increment if provided - ph.write_hydrostatic_file(meshDir) + wMax_HAMS = max(wMax, max(self.w)) - ph.write_control_file(meshDir, waterDepth=self.depth, - numFreqs=-len(self.w), minFreq=self.w[0], dFreq=np.diff(self.w[:2])[0]) + nw_HAMS = int(np.ceil(wMax_HAMS/dw_HAMS)) # ensure the upper frequency of the HAMS analysis is large enough + + ph.write_control_file(meshDir, waterDepth=self.depth, iFType=3, oFType=4, # inputs are in rad/s, outputs in s + numFreqs=-(nw_HAMS+1), minFreq=0.0, dFreq=dw_HAMS) - ph.run_hams(meshDir) # <<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< + # execute the HAMS analysis + ph.run_hams(meshDir) + # read the HAMS WAMIT-style output files data1 = os.path.join(meshDir, 'Output','Wamit_format','Buoy.1') data3 = os.path.join(meshDir, 'Output','Wamit_format','Buoy.3') - #raftDir = os.path.dirname(__file__) - #addedMass, damping = ph.read_wamit1(os.path.join(raftDir, data1)) - addedMass, damping, w_HAMS = ph.read_wamit1B(data1) - fExMod, fExPhase, fExReal, fExImag = ph.read_wamit3(data3) - - #addedMass, damping = ph.read_wamit1(data1) # original - #addedMass, damping = ph.read_wamit1('C:\\Code\\RAFT\\raft\\data\\cylinder\\Output\\Wamit_format\\Buoy.1') - #addedMass, damping = ph.read_wamit1('C:\\Code\\RAFT\\raft\\BEM\\Output\\Wamit_format\\Buoy.1') + addedMass, damping, w1 = ph.read_wamit1B(data1, TFlag=True) + M, P, R, I, w3, headings = ph.read_wamit3B(data3, TFlag=True) + + + # if requested, write hydro files for OpenFAST at whatever frequency resolution HAMS used + if len(FAST_outname) > 0: - #fExMod, fExPhase, fExReal, fExImag = ph.read_wamit3(data3) # original - #fExMod, fExPhase, fExReal, fExImag = ph.read_wamit3('C:\\Code\\RAFT\\raft\\data\\cylinder\\Output\\Wamit_format\\Buoy.3') - #fExMod, fExPhase, fExReal, fExImag = ph.read_wamit3('C:\\Code\\RAFT\\raft\\BEM\\Output\\Wamit_format\\Buoy.3') + # do a separate HAMS run to get the zero and near-infinite frequency results + ph.write_control_file(meshDir, waterDepth=self.depth, iFType=3, oFType=4, # inputs are in rad/s, outputs in s + numFreqs=-2, minFreq=0.0, dFreq=wInf) + # execute the HAMS analysis + ph.run_hams(meshDir) + + addedMassInf, dampingInf, w1Inf = ph.read_wamit1B(data1, TFlag=True) + + # radiation file + with open(FAST_outname+'.1', 'w') as f1: + + # note: in WAMIT-style output, regardless if the temporal column is period or frequency, + # a negative value indicates zero-frequency and a zero value indicates infinite frequency. + + for i in range(6): # infinite-frequency added mass (approximated from a large frequency values) + for j in range(6): + #f1.write("{:14.6e} {:d} {:d} {:13.6e}\n".format(-1.0, i+1, j+1, addedMass[i,j,-1])) + f1.write("{:14.6e} {:d} {:d} {:13.6e}\n".format(-1.0, i+1, j+1, addedMassInf[i,j,1])) + + for i in range(6): # zero-frequency added mass + for j in range(6): + f1.write("{:14.6e} {:d} {:d} {:13.6e}\n".format(0.0, i+1, j+1, addedMass[i,j,0])) + + for iw in range(1,nw_HAMS+1): + for i in range(6): # finite periods + for j in range(6): + f1.write("{:14.6e} {:d} {:d} {:13.6e} {:13.6e}\n".format( + 2.0*np.pi/w1[iw], i+1, j+1, addedMass[i,j,iw], damping[i,j,iw])) + + # excitation file + with open(FAST_outname+'.3', 'w') as f1: + + for iw in range(1,nw_HAMS+1): + for i in range(6): + for j in range(len(headings)): + f1.write("{:14.6e} {:14.6e} {:d} {:13.6e} {:13.6e} {:13.6e} {:13.6e}\n".format( + 2.0*np.pi/w1[iw], headings[j], i+1, M[j,i,iw], P[j,i,iw], R[j,i,iw], I[j,i,iw])) + + + # interpole to the frequencies RAFT is using + addedMassInterp = interp1d(w1, addedMass, assume_sorted=False, axis=2)(self.w) + dampingInterp = interp1d(w1, damping, assume_sorted=False, axis=2)(self.w) + fExRealInterp = interp1d(w3, R , assume_sorted=False )(self.w) + fExImagInterp = interp1d(w3, I , assume_sorted=False )(self.w) # copy results over to the FOWT's coefficient arrays - self.A_BEM = self.rho_water * addedMass - self.B_BEM = self.rho_water * damping - self.X_BEM = self.rho_water * self.g * (fExReal + 1j*fExImag) # linear wave excitation coefficients - self.w_BEM = w_HAMS # <<< clean this up - - # >>> do we want to seperate out infinite-frequency added mass? <<< + self.A_BEM = self.rho_water * addedMassInterp + self.B_BEM = self.rho_water * dampingInterp + self.X_BEM = self.rho_water * self.g * (fExRealInterp + 1j*fExImagInterp) + + # note: RAFT will only be using finite-frequency potential flow coefficients def calcTurbineConstants(self, case, ptfm_pitch=0): diff --git a/raft/raft_model.py b/raft/raft_model.py index bb19c0e34..08f739960 100644 --- a/raft/raft_model.py +++ b/raft/raft_model.py @@ -4,6 +4,7 @@ import numpy as np import matplotlib.pyplot as plt from matplotlib import cm +import yaml import moorpy as mp import raft.raft_fowt as fowt @@ -58,7 +59,7 @@ def __init__(self, design, nTurbines=1): self.k[i] = waveNumber(self.w[i], self.depth) # set up the FOWT here <<< only set for 1 FOWT for now <<< - self.fowtList.append(fowt.FOWT(design, w=self.w, mpb=self.ms.bodyList[0], depth=self.depth)) + self.fowtList.append(fowt.FOWT(design, self.w, self.ms.bodyList[0], depth=self.depth)) self.coords.append([0.0,0.0]) self.nDOF += 6 @@ -596,6 +597,25 @@ def calcOutputs(self): return self.results + def preprocess_HAMS(self, FAST_outname, dw=0, wMax=0): + '''This generates a mesh for the platform, runs a BEM analysis on it + using pyHAMS, and writes .1 and .3 output files for use with OpenFAST. + The mesh is only made for non-interesecting members flagged with potMod=1. + + PARAMETERS + ---------- + FAST_outname : string + The full path and file name (minus extension) for the .1 and .3 files to be written to. + dw : float + Optional specification of custom frequency increment (rad/s). + wMax : float + Optional specification of maximum frequency for BEM analysis (rad/s). Will only be + used if it is greater than the maximum frequency used in RAFT. + ''' + + self.fowtList[0].calcBEM(FAST_outname=FAST_outname, dw=dw, wMax=wMax) + + def plot(self, hideGrid=False): '''plots the whole model, including FOWTs and mooring system...''' @@ -639,7 +659,7 @@ def runRAFT(input_file, turbine_file=""): depth = float(design['mooring']['water_depth']) # for now, turn off potMod in the design dictionary to avoid BEM analysis - design['platform']['potModMaster'] = 1 + #design['platform']['potModMaster'] = 1 # read in turbine data and combine it in # if len(turbine_file) > 0: @@ -656,6 +676,8 @@ def runRAFT(input_file, turbine_file=""): model.plot() + #model.preprocess_HAMS("testHAMSoutput", dw=0.1, wMax=10) + plt.show() return model @@ -665,4 +687,5 @@ def runRAFT(input_file, turbine_file=""): import raft model = runRAFT(os.path.join(raft_dir,'designs/VolturnUS-S.yaml')) + #model = runRAFT(os.path.join(raft_dir,'designs/OC3spar.yaml')) fowt = model.fowtList[0]