minor marine growth bug fix, clean up __init__ methods

-- marine growth getMoorPyArray(pristine=False) bug fix
-- marine growth now unloads to yaml
-- clean up __init__ methods to remove unneeded inputs
-- other various small changes as recommended by Jun Tanemoto's FAModel code review (thanks Jun!!)

Author: lsirkis

famodel/cables/dynamic_cable.py

@@ -23,7 +23,7 @@ def __init__(self, id, dd=None, subsystem=None, rA=[0,0,0], rB=[0,0,0],
                  rad_anch=None, rad_fair=5, z_anch=-100, z_fair=-14, 
                  rho=1025, g=9.81,span=2000,length=2200,A=None, 
                  zJTube=-30,voltage=66,powerRating=None,cable_type=None,
-                 headingA=None,headingB=None,buoyancy_sections=None,shared=0,**kwargs):
+                 headingA=None,headingB=None,shared=0,**kwargs):
         '''
         Parameters
         ----------

famodel/cables/static_cable.py

@@ -17,11 +17,9 @@ class StaticCable(Edge):
     seabed.
     '''
 
-    def __init__(self, id, dd=None, subsystem=None, anchor=None, rA=[0,0,0], rB=[0,0,0],
-                 rad_anch=500, rad_fair=58, z_anch=-100, z_fair=-14, 
-                 rho=1025, g=9.81, L=2000,A=0,conductorSize=None, 
-                 type='static',voltage=66,powerRating=None,cable_type=None,routing=[],
-                 burial=None,zAnchor=None,**kwargs):
+    def __init__(self, id, dd=None, rA=[0,0,0], rB=[0,0,0],
+                 rho=1025, g=9.81, routing=[],
+                 burial=None, **kwargs):
         '''
         Parameters
         ----------
@@ -39,7 +37,7 @@ def __init__(self, id, dd=None, subsystem=None, anchor=None, rA=[0,0,0], rB=[0,0
         self.cableType = self.makeCableType(self.dd['cable_type'])  # Process/check it into a new dict
 
         # cable length
-        self.getLength
+        self.getLength()
 
         # end point absolute coordinates, to be set later
         self.rA = rA

famodel/helpers.py

@@ -139,14 +139,14 @@ def adjustCable(cc,project,na=None,nb=None,routeAdjustLength=500,rad_fair=None):
         hA = np.radians(90) - np.arctan2((cc.subcomponents[0].rB[0]-cc.subcomponents[0].rA[0]),(cc.subcomponents[0].rB[1]-cc.subcomponents[0].rA[1]))
     else:
         hA = np.radians(na) #headingA
-        cx[0] = [cc.attached_to[0].r[0]+500*np.cos(np.radians(na))]
-        cy[0] = [cc.attached_to[0].r[1]+500*np.sin(np.radians(na))]
+        cx[0] = [cc.attached_to[0].r[0]+routeAdjustLength*np.cos(np.radians(na))]
+        cy[0] = [cc.attached_to[0].r[1]+routeAdjustLength*np.sin(np.radians(na))]
     if nb==None:
         hB = np.radians(90) - np.arctan2((cc.subcomponents[-1].rA[0]-cc.subcomponents[-1].rB[0]),(cc.subcomponents[-1].rA[1]-cc.subcomponents[-1].rB[1]))
     else:
         hB = np.radians(nb)
-        cx[-1] = [cc.attached_to[1].r[0]+500*np.cos(np.radians(nb))]
-        cy[-1] = [cc.attached_to[1].r[1]+500*np.sin(np.radians(nb))]
+        cx[-1] = [cc.attached_to[1].r[0]+routeAdjustLength*np.cos(np.radians(nb))]
+        cy[-1] = [cc.attached_to[1].r[1]+routeAdjustLength*np.sin(np.radians(nb))]
     cc.reposition(project=project,headings=[hA,hB],rad_fair=rad_fair)
 
 # find cable(s) associated with specific platform(s) xy coordinates
@@ -251,6 +251,8 @@ def head_adjust(att,heading,rad_buff=np.radians(30),endA_dir=1, adj_dir=1):
     '''
     if heading<0:
         headnew = np.pi*2 + heading
+    elif heading>2*np.pi:
+        heading - 2*np.pi
     else:
         headnew = heading
     attheadings = [] # complete list of mooring headings to avoid, from all platforms

famodel/mooring/mooring.py

@@ -806,7 +806,7 @@ def ddSave(self):
         return self.dd
     """
 
-    def addMarineGrowth(self, mgDict, project=None, idx=None):
+    def addMarineGrowth(self, mgDict):
         '''Re-creates sections part of design dictionary to account for marine 
         growth on the subystem, then calls createSubsystem() to recreate the line
 
@@ -823,17 +823,6 @@ def addMarineGrowth(self, mgDict, project=None, idx=None):
                                                 - [0.20,-40,0]
                 rho : list of densities for each thickness, or one density for all thicknesses, [kg/m^3] (optional - default is 1325 kg/m^3)
                 }
-        project : object, optional
-            A FAModel project object, with the mooringListPristine used as the basis
-            to build the marine growth model, necessary if the addMarineGrowth method
-            will be called in a loop (or even just multiple times) to improve accuracy 
-            of change depths, which may decrease significantly after solveEquilibrium() 
-            for the moorpy model. The default is None.
-        idx : tuple, optional
-            A key for the pristineMooringList in the project object that is associated
-            with this mooring object. Since the pristineMooringList is a deepcopy of the 
-            project mooringList, the mooring objects are not the same and therefore if the 
-            project object is provided in the method call, the index must also be provided.
 
         Returns
         -------

famodel/project.py

@@ -97,6 +97,7 @@ def __init__(self, lon=0, lat=0, file=None, depth=202,raft=1):
         self.rho_water = 1025 # density of water (default to saltwater) [kg/m^3]
         self.rho_air = 1.225 # density of air [kg/m^3]
         self.mu_air = 1.81e-5 # dynamic viscosity of air [Pa*s]
+        self.marine_growth = None
 
         # Project boundary (vertical stack of x,y coordinate pairs [m])
         self.boundary = np.zeros([0,2])
@@ -110,10 +111,9 @@ def __init__(self, lon=0, lat=0, file=None, depth=202,raft=1):
         self.depth = depth
         self.grid_depth  = np.array([[self.depth]])  # depth at each grid point [iy, ix]
 
-        self.seabed_type = 'clay'  # switch of which soil property set to use ('clay', 'sand', or 'rock')
-        
         # soil parameters at each grid point
         self.soilProps = {}
+        self.soil_names = []
         self.soil_mode  = 0                # soil/anchor model level to use (0: none; 1: simple categories; 2: quantitative)
         self.soil_class = [["none"]]       # soil classification name ('clay', 'sand', or 'rock' with optional modifiers)
         self.soil_gamma = np.zeros((1,1))  # soil effective unit weight [kPa] (all soils)
@@ -1858,18 +1858,19 @@ def addAnchor(self, id=None, dd=None, design=None, cost=0, atype=None, platform=
         return(anchor)
 
     def addTurbine(self,id=None, typeID=None, platform=None, turbine_dd={}):
+        
         if typeID == None:
             typeID = len(self.turbineTypes)
         if id==None:
             id = 'T'+str(typeID)+'_'+str(len(self.turbineList))
 
         if turbine_dd and 'blade' in turbine_dd:
-            blade_diameter = turbine_dd['blade']['Rtip']*2        
+            rotor_diameter = turbine_dd['blade']['Rtip']*2        
             self.turbineTypes.append(turbine_dd)
         else:
             blade_diameter = 0
 
-        self.turbineList[id] = Turbine(turbine_dd,id,D=blade_diameter)
+        self.turbineList[id] = Turbine(turbine_dd,id,D=rotor_diameter)
         self.turbineList[id].dd['type'] = typeID
         if platform != None:
             platform.attach(self.turbineList[id])
@@ -2184,7 +2185,7 @@ def plot2d(self, ax=None, plot_seabed=False, plot_soil=False,
         depth_vmin = kwargs.get('depth_vmin', None)
         depth_vmax = kwargs.get('depth_vmax', None)
         bath_levels = kwargs.get('bath_levels', None)
-        show_legend = kwargs.get('show_legend', True)
+        plot_legend = kwargs.get('plot_legend', True)
         legend_x = kwargs.get('legend_x', 0.5)
         legend_y = kwargs.get('legend_y', -0.1)
 
@@ -2302,16 +2303,21 @@ def plot2d(self, ax=None, plot_seabed=False, plot_soil=False,
                                 line.type['material'][1:]+' Mooring')
 
                 if mooring.ss:
-                    mooring.ss.drawLine2d(0, ax, color="self", endpoints=False, 
-                                          Xuvec=[1,0,0], Yuvec=[0,1,0], line_depth_settings=line_depth_settings, label=labs)  
+                    mooring.ss.drawLine2d(0, ax, color="self", 
+                                          plot_endpoints=False, 
+                                          Xuvec=[1,0,0], Yuvec=[0,1,0], 
+                                          line_depth_settings=line_depth_settings, 
+                                          label=labs)  
                 elif mooring.parallels:
                     for i,line in enumerate(lineList):
                         line.drawLine2d(0, ax, color="self",
-                                        Xuvec=[1,0,0], Yuvec=[0,1,0], label=labs[i])
+                                        Xuvec=[1,0,0], Yuvec=[0,1,0], 
+                                        label=labs[i])
 
                 else: # simple line plot
                     ax.plot([mooring.rA[0], mooring.rB[0]], 
-                            [mooring.rA[1], mooring.rB[1]], 'k', lw=0.5, label='Mooring Line')
+                            [mooring.rA[1], mooring.rB[1]], 'k', lw=0.5, 
+                            label='Mooring Line')
 
         # ---- Add colorbar for line depth ----
         if line_depth_settings is not None and not bare:
@@ -2351,23 +2357,23 @@ def plot2d(self, ax=None, plot_seabed=False, plot_soil=False,
                             # has routing  - first plot rA to sub.coordinate[0] connection
                             ax.plot([sub.rA[0],sub.x[0]],
                                     [sub.rA[1],sub.y[0]],':',color = Ccable,
-                                    lw=1.2,label='Static Cable '+str(cableSize)+' mm$^{2}$')
+                                    lw=1.2,label=f'Static Cable {cableSize} mm$^{2}$')
                             # now plot route
                             if len(sub.x) > 1:
                                 for i in range(1,len(sub.x)):
                                     ax.plot([sub.x[i-1],sub.x[i]],
                                             [sub.y[i-1],sub.y[i]],
                                             ':', color=Ccable, lw=1.2,
-                                            label='Static Cable '+str(cableSize)+' mm$^{2}$')
+                                            label=f'Static Cable {cableSize} mm$^{2}$')
                             # finally plot sub.coordinates[-1] to rB connection
                             ax.plot([sub.x[-1],sub.rB[0]],
                                     [sub.y[-1],sub.rB[1]],':',color=Ccable,
-                                    lw=1.2,label='Static Cable '+str(cableSize)+' mm$^{2}$')
+                                    lw=1.2,label=f'Static Cable {cableSize} mm$^{2}$')
                         else:
                             # if not routing just do simple line plot
                             ax.plot([sub.rA[0],sub.rB[0]], 
                                     [sub.rA[1], sub.rB[1]],':',color = Ccable, lw=1.2,
-                                    label='Static Cable '+str(cableSize)+' mm$^{2}$')
+                                    label=f'Static Cable {cableSize} mm$^{2}$')
 
                         # if cable_labels:
                         #     x = np.mean([sub.rA[0],sub.rB[0]])
@@ -2380,7 +2386,7 @@ def plot2d(self, ax=None, plot_seabed=False, plot_soil=False,
                     elif isinstance(sub,DynamicCable):
                             ax.plot([sub.rA[0],sub.rB[0]], 
                                     [sub.rA[1], sub.rB[1]],color = Ccable, lw=1.2,
-                                    label='Dynamic Cable '+str(cableSize)+' mm$^{2}$')
+                                    label=f'Dynamic Cable {cableSize} mm$^{2}$')
 
                             if cable_labels:
                                 x = np.mean([sub.rA[0],sub.rB[0]])
@@ -2421,10 +2427,17 @@ def plot2d(self, ax=None, plot_seabed=False, plot_soil=False,
             labels += [h.get_label() for h in soil_handles]
         by_label = dict(zip(labels, handles))  # Removing duplicate labels
 
-        if show_legend:
+        if plot_legend:
             ax.legend(by_label.values(), by_label.keys(),loc='upper center',bbox_to_anchor=(legend_x, legend_y), fancybox=True, ncol=4)
         if save:
-            plt.savefig('2dfarm.png', dpi=600, bbox_inches='tight')  # Adjust the dpi as needed
+            counter = 1
+            output_filename = f'2dfarm_{counter}.png'
+            while os.path.exists(output_filename):
+                counter += 1
+                output_filename = f'2dfarm_{counter}.png'
+            
+            # Increase the resolution when saving the plot
+            plt.savefig(output_filename, dpi=600, bbox_inches='tight')  # Adjust the dpi as needed
 
             # TODO - add ability to plot from RAFT FOWT
         if return_contour:
@@ -2483,21 +2496,8 @@ def plot3d(self, ax=None, figsize=(10,8), plot_fowt=False, save=False,
                     args_bath = {'cmap': cmap, 'vmin':min([min(x) for x in -self.grid_depth]),
                                  'vmax': vmax}
 
-                # >>> below code sections seem to do unnecessary interpolation and regridding... >>>
                 if boundary_only:   # if you only want to plot the bathymetry that's underneath the boundary, rather than the whole file
-                    boundary = np.vstack([self.boundary, self.boundary[0,:]])
-                    xs = np.linspace(min(boundary[:,0]),max(boundary[:,0]),len(boundary[:,0]))
-                    ys = np.linspace(min(boundary[:,1]),max(boundary[:,1]),len(boundary[:,1]))
-
-                    self.setGrid(xs, ys)
-            
-                    # plot the bathymetry in matplotlib using a plot_surface
-                    X, Y = np.meshgrid(xs, ys)  # 2D mesh of seabed grid
-
-                    plot_depths = np.zeros([len(ys), len(xs)])
-                    for i in range(len(ys)):
-                        for j in range(len(xs)):
-                            plot_depths[i,j], nvec = self.getDepthAtLocation(xs[j], ys[i])
+                    self.trimGrids()
 
                 else:
                     xs = np.linspace(min(self.grid_x),max(self.grid_x),len(self.grid_x))
@@ -2775,9 +2775,12 @@ def getMoorPyArray(self, plt=0, pristineLines=True, cables=True):
                     # set location of subsystem for simpler coding
                     ssloc = mooring.ss
                 else:
-                    mooring.createSubsystem(pristine=False, ms=self.ms)
+                    # if modified , assume the ss_mod is already generated
                     # set location of subsystem for simpler coding
                     ssloc = mooring.ss_mod
+                    # assign ss_mod in to ms
+                    self.ms.lineList.append(ssloc)
+                    ssloc.number = len(self.ms.lineList)
 
                 # (ms.lineList.append is now done in Mooring.createSubsystem)
 
@@ -2889,19 +2892,19 @@ def getMoorPyArray(self, plt=0, pristineLines=True, cables=True):
             if check == 1: # mooring object not in any anchor lists
                 # new shared line
                 # create subsystem for shared line
-                if hasattr(mooring, 'shared'):  # <<<
-                    mooring.createSubsystem(case=mooring.shared, 
-                        pristine=pristineLines, ms=self.ms)
-                else:
-                    mooring.createSubsystem(case=1,pristine=pristineLines, 
-                        ms=self.ms) # we doubled all symmetric lines so any shared lines should be case 1
+   
                 # set location of subsystem for simpler coding
                 if pristineLines:
+                    mooring.createSubsystem(case=mooring.shared, 
+                        pristine=pristineLines, ms=self.ms)
                     ssloc = mooring.ss
                 else:
+                    # modified line should already exist - just attach it
                     ssloc = mooring.ss_mod
+                    self.ms.lineList.append(ssloc)
+                    ssloc.number = len(self.ms.lineList)
 
-                # (ms.lineList.append is now done in Mooring.createSubsystem)
+                # (ms.lineList.append is now done in Mooring.createSubsystem unless not pristine)
 
                 # find associated platforms/ buoys
                 att = mooring.attached_to
@@ -3432,7 +3435,7 @@ def getMarineGrowth(self,mgDict_start=None,lines='all',tol=2,display=False):
                     for j in range(0,len(cP)):
                         cEq.append(mgDict_start['th'][cD[j][0]][cD[j][1]] - self.cableList[i[0]].subcomponents[i[1]].ss_mod.pointList[cP[j]].r[2])
                 else:
-                    cD,cP = self.mooringList[i].addMarineGrowth(mgDict,project=self,idx=i)
+                    cD,cP = self.mooringList[i].addMarineGrowth(mgDict)
                     for j in range(0,len(cP)):
                         cEq.append(mgDict_start['th'][cD[j][0]][cD[j][1]] - self.mooringList[i].ss_mod.pointList[cP[j]].r[2])
                 # adjust depth to change based on difference between actual and desired change depth
@@ -4532,6 +4535,14 @@ def unload(self,file='project.yaml'):
             site['bathymetry'] = {'x':[float(x) for x in self.grid_x],'y':[float(y) for y in self.grid_y],'depths':[[float(y) for y in x] for x in self.grid_depth]}
         if len(self.boundary)>0:
             site['boundaries'] = {'x_y':[[float(y) for y in x] for x in self.boundary]}
+            
+        if self.marine_growth:
+            site['marine_growth'] = {
+                'keys':['thickness', 'lowerRange', 'upperRange', 'density'][:len(self.marine_growth['th'])],
+                'data':[v for v in self.marine_growth['th']]
+                }
+            if 'buoy_th' in self.marine_growth:
+                site['marine_growth']['buoys'] = [x for x in self.marine_growth['buoy_th']]
         site['general'] = {'water_depth':float(self.depth),'rho_air':float(self.rho_air),
                            'rho_water':float(self.rho_water),'mu_air':float(self.mu_air)}
 
@@ -4870,7 +4881,10 @@ def extractFarmInfo(self, cmax=5, fmax=10/6, Cmeander=1.9):
                 phi_deg       = (phi_deg + 180) % 360 - 180  # Shift range to -180 to 180
                 for att in pf.attachments.values():
                     if isinstance(att['obj'],Turbine):
-                        D    = 240   # att['obj'].D         (assuming 15MW)
+                        if hasattr(att['obj'], 'D'):
+                            D = att['obj'].D
+                        else:
+                            D = 242
                         zhub = att['obj'].dd['hHub']
 
                 wts[i] = {