StFile

wetb/hawc2/log_file.py

@@ -168,7 +168,7 @@ class LogFile(LogInterpreter):
                 fid.seek(self.position)
                 txt = fid.read()
             self.position += len(txt)
-            txt = txt.decode(encoding='utf_8', errors='strict')
+            txt = txt.decode(encoding='cp1252', errors='strict')
             if txt != "":
                 LogInterpreter.update_status(self, txt)
 

wetb/hawc2/simulation.py

@@ -15,7 +15,7 @@ import subprocess
 import sys
 import time
 from wetb.hawc2 import log_file
-from wetb.hawc2.htc_file import HTCFile
+from wetb.hawc2.htc_file import HTCFile, H2aeroHTCFile
 from wetb.hawc2.log_file import LogFile, LogInfo
 
 from future import standard_library
@@ -26,6 +26,7 @@ from wetb.utils.cluster_tools.cluster_resource import LocalResource
 from wetb.utils.cluster_tools.pbsjob import SSHPBSJob, DONE, NOT_SUBMITTED
 from wetb.utils.cluster_tools.ssh_client import SSHClient
 from wetb.utils.timing import print_time
+from _datetime import datetime
 standard_library.install_aliases()
 from threading import Thread
 
@@ -84,7 +85,7 @@ class Simulation(object):
         if not os.path.isabs(htcfilename):
             htcfilename = os.path.join(modelpath, htcfilename)
         self.filename = os.path.basename(htcfilename)
-        self.htcFile = HTCFile(htcfilename)
+        self.htcFile = H2aeroHTCFile(htcfilename)
         self.time_stop = self.htcFile.simulation.time_stop[0]
         self.hawc2exe = hawc2exe
         self.copy_turbulence = copy_turbulence
@@ -174,7 +175,7 @@ class Simulation(object):
             assert not src.startswith(".."), "%s referes to a file outside the model path\nAll input files be inside model path" % src
             return src
         input_patterns = [fmt(src) for src in self.htcFile.input_files() + self.htcFile.turbulence_files() + self.additional_files().get('input', [])]
-        input_files = set([f for pattern in input_patterns for f in glob.glob(os.path.join(self.modelpath, pattern))])
+        input_files = set([f for pattern in input_patterns for f in glob.glob(os.path.join(self.modelpath, pattern)) if os.path.isfile(f)])
         if not os.path.isdir(os.path.dirname(self.modelpath + self.stdout_filename)):
             os.makedirs(os.path.dirname(self.modelpath + self.stdout_filename))
         self.host._prepare_simulation(input_files)
@@ -190,6 +191,7 @@ class Simulation(object):
 
     def simulate(self):
         #starts blocking simulation
+        #self.simulation_start_time = self.host.get_datetime()
         self.is_simulating = True
         self.errors = []
         self.status = INITIALIZING
@@ -376,6 +378,8 @@ class LocalSimulationHost(SimulationResource):
         self.resource = resource
         self.simulationThread = SimulationThread(self.sim)
 
+    def get_datetime(self):
+        return datetime.now()
 
     def glob(self, path):
         return glob.glob(path)
@@ -478,7 +482,7 @@ class SimulationThread(Thread):
         self.process.communicate()
         errorcode = self.process.returncode
 
-        with open(self.modelpath + self.sim.stdout_filename, encoding='utf-8') as fid:
+        with open(self.modelpath + self.sim.stdout_filename, encoding='cp1252') as fid:
             stdout = fid.read()
         self.res = errorcode, stdout
 
@@ -497,6 +501,11 @@ class PBSClusterSimulationHost(SimulationResource, SSHClient):
     hawc2exe = property(lambda self : os.path.basename(self.sim.hawc2exe))
 
 
+    def get_datetime(self):
+        v, out, err = self.execute('date "+%Y,%m,%d,%H,%M,%S"')
+        if v == 0:
+            return datetime.strptime(out.strip(), "%Y,%m,%d,%H,%M,%S")
+
     def _prepare_simulation(self, input_files):
         with self:
             self.execute(["mkdir -p .hawc2launcher/%s" % self.simulation_id], verbose=False)
@@ -517,6 +526,7 @@ class PBSClusterSimulationHost(SimulationResource, SSHClient):
 
 
 
+
     def _finish_simulation(self, output_files):
         with self:
             for src_file in output_files:

wetb/hawc2/st_file.py

+'''
+Created on 24/04/2014
+
+@author: MMPE
+'''
+from __future__ import print_function
+from __future__ import unicode_literals
+from __future__ import division
+from __future__ import absolute_import
+from io import open
+from builtins import range
+from builtins import int
+from future import standard_library
+standard_library.install_aliases()
+
+import numpy as np
+
+class StFile(object):
+    """Read HAWC2 St (beam element structural data) file
+
+    Methods are autogenerated for:
+
+    - r : curved length distance from main_body node 1 [m]
+    - m : mass per unit length [kg/m]
+    - xm : xc2-coordinate from C1/2 to mass center [m]
+    - ym : yc2-coordinate from C1/2 to mass center [m]
+    - rix : radius of gyration related to elastic center. Corresponds to rotation about principal bending xe axis [m]
+    - riy : radius of gyration related to elastic center. Corresponds to rotation about principal bending ye axis [m]
+    - xs : xc2-coordinate from C1/2 to shear center [m]. The shear center is the point where external forces only contributes to pure bending and no torsion.
+    - ys : yc2-coordinate from C1/2 to shear center [m]. The shear center is the point where external forces only contributes to pure bending and no torsion.
+    - E : modulus of elasticity [N/m2]
+    - G : shear modulus of elasticity [N/m2]
+    - Ix : area moment of inertia with respect to principal bending xe axis [m4]. This is the principal bending axis most parallel to the xc2 axis
+    - Iy : area moment of inertia with respect to principal bending ye axis [m4]
+    - K : torsional stiffness constant with respect to ze axis at the shear center [m4/rad]. For a circular section only this is identical to the polar moment of inertia.
+    - kx : shear factor for force in principal bending xe direction [-]
+    - ky : shear factor for force in principal bending ye direction [-]
+    - A : cross sectional area [m2]
+    - pitch : structural pitch about z_c2 axis. This is the angle between the xc2 -axis defined with the c2_def command and the main principal bending axis xe.
+    - xe : xc2-coordinate from C1/2 to center of elasticity [m]. The elastic center is the point where radial force (in the z-direction) does not contribute to bending around the x or y directions.
+    - ye : yc2-coordinate from C1/2 to center of elasticity [m]. The elastic center is the point where radial force (in the
+
+    The autogenerated methods have the following structure
+
+    def xxx(radius=None, mset_nr=1, set_nr=1):
+        Parameters:
+        -----------
+        radius : int, float, array_like or None, optional
+            Radius/radii of interest\n
+            If int, float or array_like: values are interpolated to requested radius/radii
+            If None (default): Values of all radii specified in st file returned
+        mset_nr : int, optional
+            Main set number
+        set_nr : int, optional
+            Sub set number
+
+
+    Examples
+    --------
+    >>> stfile = StFile(r"tests/test_files/DTU_10MW_RWT_Blade_st.dat")
+    >>> print (stfile.m()) # Interpolated mass at radius 36
+    [ 1189.51054664  1191.64291781  1202.76694262  ... 15.42438683]
+    >>> print (st.E(radius=36, mset_nr=1, set_nr=1))  # Elasticity interpolated to radius 36m
+    8722924514.652649
+    >>> print (st.E(radius=36, mset_nr=1, set_nr=2))  # Same for stiff blade set
+    8.722924514652648e+17
+    """
+    def __init__(self, filename):
+        with open (filename) as fid:
+            txt = fid.read()
+        no_maindata_sets = int(txt.strip()[0])
+        assert no_maindata_sets == txt.count("#")
+        self.main_data_sets = {}
+        for mset in txt.split("#")[1:]:
+            mset_nr = int(mset.strip().split()[0])
+            set_data_dict = {}
+
+            for set_txt in mset.split("$")[1:]:
+                set_lines = set_txt.split("\n")
+                set_nr, no_rows = map(int, set_lines[0].split()[:2])
+                assert set_nr not in set_data_dict
+                set_data_dict[set_nr] = np.array([set_lines[i].split() for i in range(1, no_rows + 1)], dtype=np.float)
+            self.main_data_sets[mset_nr] = set_data_dict
+
+        for i, name in enumerate("r m x_cg y_cg ri_x ri_y x_sh y_sh E G I_x I_y I_p k_x k_y A pitch x_e y_e".split()):
+            setattr(self, name, lambda radius=None, mset_nr=1, set_nr=1, column=i: self._value(radius, column, mset_nr, set_nr))
+
+    def _value(self, radius, column, mset_nr=1, set_nr=1):
+        st_data = self.main_data_sets[mset_nr][set_nr]
+        if radius is None:
+            radius = self.radius(None, mset_nr, set_nr)
+        return np.interp(radius, st_data[:, 0], st_data[:, column])
+
+    def radius(self, radius=None, mset_nr=1, set_nr=1):
+        r = self.main_data_sets[mset_nr][set_nr][:, 0]
+        if radius is None:
+            return r
+        return r[np.argmin(np.abs(r - radius))]
+
+
+
+if __name__ == "__main__":
+    st = StFile(r"C:\mmpe\HAWC2\models\DTU10MWRef6.0\data/DTU_10MW_RWT_Blade_st.dat")
+    print (st.m())
+    print (st.E(radius=36, mset_nr=1, set_nr=1))  # Elastic blade
+    print (st.E(radius=36, mset_nr=1, set_nr=2))  # stiff blade
+    #print (st.radius())
+    xyz = np.array([st.x_e(), st.y_e(), st.r()]).T[:40]
+    n = 2
+    xyz = np.array([st.x_e(None, 1, n), st.y_e(None, 1, n), st.r(None, 1, n)]).T[:40]
+    #print (xyz)
+    print (np.sqrt(np.sum((xyz[1:] - xyz[:-1]) ** 2, 1)).sum())
+    print (xyz[-1, 2])
+    print (np.sqrt(np.sum((xyz[1:] - xyz[:-1]) ** 2, 1)).sum() - xyz[-1, 2])
+    print (st.x_e(67.8883), st.y_e(67.8883))
+    #print (np.sqrt(np.sum(np.diff(xyz, 0) ** 2, 1)))
+    print (st.pitch(67.8883 - 0.01687))
+    print (st.pitch(23.2446))
+
+
+
+    #print (st.)
+    #print (st.)