<s> import copy <EOL> import numpy as np <EOL> from bspline import Bspline <EOL> class Coordinates ( object ) : <EOL> """<STR_LIT>""" <EOL> def __init__ ( self , points , cartesian = True ) : <EOL> """<STR_LIT>""" <EOL> points = np . array ( points , dtype = np . float64 ) <EOL> if cartesian : <EOL> self . cartesian = points <EOL> X , Y , Z = points [ : , <NUM_LIT:0> ] , points [ : , <NUM_LIT:1> ] , points [ : , <NUM_LIT:2> ] <EOL> R = np . sqrt ( Y ** <NUM_LIT:2> + Z ** <NUM_LIT:2> ) <EOL> Theta = np . arctan2 ( Y , Z ) <EOL> self . cylindrical = np . nan_to_num ( np . vstack ( ( X , R , Theta ) ) . T ) <EOL> else : <EOL> self . cylindrical = points <EOL> X , R , Theta = points [ : , <NUM_LIT:0> ] , points [ : , <NUM_LIT:1> ] , points [ : , <NUM_LIT:2> ] <EOL> Z = R * np . cos ( Theta ) <EOL> Y = R * np . sin ( Theta ) <EOL> self . cartesian = np . nan_to_num ( np . vstack ( ( X , Y , Z ) ) . T ) <EOL> class Body ( object ) : <EOL> """<STR_LIT>""" <EOL> def __init__ ( self , stl , controls , name = "<STR_LIT:body>" , r_ref = None , x_ref = None ) : <EOL> """<STR_LIT>""" <EOL> self . stl = stl <EOL> geom_points = stl . points <EOL> self . coords = Coordinates ( geom_points , cartesian = True ) <EOL> self . P = self . coords . cylindrical <EOL> self . P_cart = self . coords . cartesian <EOL> self . P_bar = geom_points . copy ( ) <EOL> if isinstance ( controls , int ) : <EOL> X = geom_points [ : , <NUM_LIT:0> ] <EOL> x_max = np . max ( X ) <EOL> x_min = np . min ( X ) <EOL> C_x = np . linspace ( x_min , x_max , controls ) <EOL> C_r = np . zeros ( ( controls , ) ) <EOL> control_points = np . array ( zip ( C_x , C_r ) ) <EOL> self . C = control_points <EOL> self . n_controls = controls <EOL> else : <EOL> self . C = controls <EOL> self . n_controls = len ( control_points ) <EOL> self . C_bar = self . C . copy ( ) <EOL> self . delta_C = np . zeros ( self . C . shape ) <EOL> self . bs = Bspline ( self . C , geom_points ) <EOL> self . name = name <EOL> if x_ref is not None : <EOL> self . x_mag = float ( x_ref ) <EOL> else : <EOL> self . x_mag = <NUM_LIT:10> ** np . floor ( np . log10 ( np . average ( geom_points [ : , <NUM_LIT:0> ] ) ) ) <EOL> if r_ref is not None : <EOL> self . r_mag = float ( r_ref ) <EOL> else : <EOL> indecies = np . logical_and ( abs ( geom_points [ : , <NUM_LIT:2> ] ) < <NUM_LIT> , geom_points [ : , <NUM_LIT:1> ] > <NUM_LIT:0> ) <EOL> points = geom_points [ indecies ] <EOL> self . r_mag = <NUM_LIT:10> ** np . floor ( np . log10 ( np . average ( points [ : , <NUM_LIT:1> ] ) ) ) <EOL> self . Theta = self . P [ : , <NUM_LIT:2> ] <EOL> self . sin_theta = np . tile ( np . sin ( self . Theta ) , ( self . n_controls , <NUM_LIT:1> ) ) . T . flatten ( ) <EOL> self . cos_theta = np . tile ( np . cos ( self . Theta ) , ( self . n_controls , <NUM_LIT:1> ) ) . T . flatten ( ) <EOL> self . dP_bar_xqdC = np . array ( self . x_mag * self . bs . B . flatten ( ) ) <EOL> self . dP_bar_rqdC = np . array ( self . r_mag * self . bs . B . flatten ( ) ) <EOL> self . dXqdC = self . dP_bar_xqdC . reshape ( - <NUM_LIT:1> , self . n_controls ) <EOL> self . dYqdC = ( self . dP_bar_rqdC * self . sin_theta ) . reshape ( - <NUM_LIT:1> , self . n_controls ) <EOL> self . dZqdC = ( self . dP_bar_rqdC * self . cos_theta ) . reshape ( - <NUM_LIT:1> , self . n_controls ) <EOL> def copy ( self ) : <EOL> return copy . deepcopy ( self ) <EOL> def deform ( self , delta_C ) : <EOL> """<STR_LIT>""" <EOL> self . delta_C = delta_C <EOL> self . delta_C [ : , <NUM_LIT:0> ] = self . delta_C [ : , <NUM_LIT:0> ] * self . x_mag <EOL> self . C_bar = self . C + self . delta_C <EOL> delta_P = self . bs . calc ( self . C_bar ) <EOL> self . P_bar = self . P . copy ( ) <EOL> self . P_bar [ : , <NUM_LIT:0> ] = delta_P [ : , <NUM_LIT:0> ] <EOL> self . P_bar [ : , <NUM_LIT:1> ] = self . P [ : , <NUM_LIT:1> ] + self . r_mag * delta_P [ : , <NUM_LIT:1> ] <EOL> self . coords = Coordinates ( self . P_bar , cartesian = False ) <EOL> self . P_bar_cart = self . coords . cartesian <EOL> self . Xo = self . P_bar_cart [ : , <NUM_LIT:0> ] <EOL> self . Yo = self . P_bar_cart [ : , <NUM_LIT:1> ] <EOL> self . Zo = self . P_bar_cart [ : , <NUM_LIT:2> ] <EOL> self . stl . update_points ( self . P_bar_cart ) <EOL> return self . P_bar <EOL> class Shell ( object ) : <EOL> """<STR_LIT>""" <EOL> def __init__ ( self , outer_stl , inner_stl , center_line_controls , <EOL> thickness_controls , name = '<STR_LIT>' , r_ref = None , x_ref = None ) : <EOL> self . outer_stl = outer_stl <EOL> self . inner_stl = inner_stl <EOL> outer_points = outer_stl . points <EOL> inner_points = inner_stl . points <EOL> self . n_outer = len ( outer_points ) <EOL> self . n_inner = len ( inner_points ) <EOL> self . outer_coords = Coordinates ( outer_points , cartesian = True ) <EOL> self . inner_coords = Coordinates ( inner_points , cartesian = True ) <EOL> self . Po = self . outer_coords . cylindrical <EOL> self . Pi = self . inner_coords . cylindrical <EOL> self . Po_cart = self . outer_coords . cartesian <EOL> self . Pi_cart = self . inner_coords . cartesian <EOL> self . Po_bar = outer_points . copy ( ) <EOL> self . Pi_bar = inner_points . copy ( ) <EOL> self . name = name <EOL> if isinstance ( center_line_controls , int ) : <EOL> X = outer_points [ : , <NUM_LIT:0> ] <EOL> x_max = np . max ( X ) <EOL> x_min = np . min ( X ) <EOL> C_x = np . linspace ( x_min , x_max , center_line_controls ) <EOL> C_r = np . zeros ( ( center_line_controls , ) ) <EOL> control_points = np . array ( zip ( C_x , C_r ) ) <EOL> self . Cc = control_points <EOL> self . n_c_controls = center_line_controls <EOL> else : <EOL> self . Cc = center_line_controls <EOL> self . n_c_controls = len ( center_line_controls ) <EOL> self . Cc_bar = self . Cc . copy ( ) <EOL> self . delta_Cc = np . zeros ( self . Cc . shape ) <EOL> if isinstance ( thickness_controls , int ) : <EOL> X = inner_points [ : , <NUM_LIT:0> ] <EOL> x_max = np . max ( X ) <EOL> x_min = np . min ( X ) <EOL> C_x = np . linspace ( x_min , x_max , thickness_controls ) <EOL> C_r = np . zeros ( ( thickness_controls , ) ) <EOL> control_points = np . array ( zip ( C_x , C_r ) ) <EOL> self . Ct = control_points <EOL> self . n_t_controls = thickness_controls <EOL> else : <EOL> self . Ct = thickness_controls <EOL> self . n_t_controls = len ( thickness_controls ) <EOL> self . Ct_bar = self . Ct . copy ( ) <EOL> self . delta_Ct = np . zeros ( self . Ct . shape ) <EOL> self . bsc_o = Bspline ( self . Cc , outer_points ) <EOL> self . bsc_i = Bspline ( self . Cc , inner_points ) <EOL> self . bst_o = Bspline ( self . Ct , outer_points ) <EOL> self . bst_i = Bspline ( self . Ct , inner_points ) <EOL> self . name = name <EOL> if x_ref is not None : <EOL> self . x_mag = float ( x_ref ) <EOL> else : <EOL> self . x_mag = <NUM_LIT:10> ** np . floor ( np . log10 ( np . average ( outer_points [ : , <NUM_LIT:0> ] ) ) ) <EOL> if r_ref is not None : <EOL> self . r_mag = float ( r_ref ) <EOL> else : <EOL> indecies = np . logical_and ( abs ( outer_points [ : , <NUM_LIT:2> ] ) < <NUM_LIT> , outer_points [ : , <NUM_LIT:1> ] > <NUM_LIT:0> ) <EOL> points = outer_points [ indecies ] <EOL> self . r_mag = <NUM_LIT:10> ** np . floor ( np . log10 ( np . average ( points [ : , <NUM_LIT:1> ] ) ) ) <EOL> self . outer_theta = self . Po [ : , <NUM_LIT:2> ] <EOL> self . sin_outer_c_theta = np . tile ( np . sin ( self . outer_theta ) , ( self . n_c_controls , <NUM_LIT:1> ) ) . T . flatten ( ) <EOL> self . cos_outer_c_theta = np . tile ( np . cos ( self . outer_theta ) , ( self . n_c_controls , <NUM_LIT:1> ) ) . T . flatten ( ) <EOL> self . sin_outer_t_theta = np . tile ( np . sin ( self . outer_theta ) , ( self . n_t_controls , <NUM_LIT:1> ) ) . T . flatten ( ) <EOL> self . cos_outer_t_theta = np . tile ( np . cos ( self . outer_theta ) , ( self . n_t_controls , <NUM_LIT:1> ) ) . T . flatten ( ) <EOL> self . inner_theta = self . Pi [ : , <NUM_LIT:2> ] <EOL> self . sin_inner_c_theta = np . tile ( np . sin ( self . inner_theta ) , ( self . n_c_controls , <NUM_LIT:1> ) ) . T . flatten ( ) <EOL> self . cos_inner_c_theta = np . tile ( np . cos ( self . inner_theta ) , ( self . n_c_controls , <NUM_LIT:1> ) ) . T . flatten ( ) <EOL> self . sin_inner_t_theta = np . tile ( np . sin ( self . inner_theta ) , ( self . n_t_controls , <NUM_LIT:1> ) ) . T . flatten ( ) <EOL> self . cos_inner_t_theta = np . tile ( np . cos ( self . inner_theta ) , ( self . n_t_controls , <NUM_LIT:1> ) ) . T . flatten ( ) <EOL> self . dPo_bar_xqdCc = np . array ( self . x_mag * self . bsc_o . B . flatten ( ) ) <EOL> self . dPo_bar_rqdCc = np . array ( self . r_mag * self . bsc_o . B . flatten ( ) ) <EOL> self . dPi_bar_xqdCc = np . array ( self . x_mag * self . bsc_i . B . flatten ( ) ) <EOL> self . dPi_bar_rqdCc = np . array ( self . r_mag * self . bsc_i . B . flatten ( ) ) <EOL> self . dPo_bar_rqdCt = np . array ( self . r_mag * self . bst_o . B . flatten ( ) ) <EOL> self . dPi_bar_rqdCt = - <NUM_LIT:1> * np . array ( self . r_mag * self . bst_i . B . flatten ( ) ) <EOL> self . dXoqdCc = self . dPo_bar_xqdCc . reshape ( - <NUM_LIT:1> , self . n_c_controls ) <EOL> self . dYoqdCc = ( self . dPo_bar_rqdCc * self . sin_outer_c_theta ) . reshape ( - <NUM_LIT:1> , self . n_c_controls ) <EOL> self . dZoqdCc = ( self . dPo_bar_rqdCc * self . cos_outer_c_theta ) . reshape ( - <NUM_LIT:1> , self . n_c_controls ) <EOL> self . dXiqdCc = self . dPi_bar_xqdCc . reshape ( - <NUM_LIT:1> , self . n_c_controls ) <EOL> self . dYiqdCc = ( self . dPi_bar_rqdCc * self . sin_inner_c_theta ) . reshape ( - <NUM_LIT:1> , self . n_c_controls ) <EOL> self . dZiqdCc = ( self . dPi_bar_rqdCc * self . cos_inner_c_theta ) . reshape ( - <NUM_LIT:1> , self . n_c_controls ) <EOL> self . dYoqdCt = ( self . dPo_bar_rqdCt * self . sin_outer_t_theta ) . reshape ( - <NUM_LIT:1> , self . n_t_controls ) <EOL> self . dZoqdCt = ( self . dPo_bar_rqdCt * self . cos_outer_t_theta ) . reshape ( - <NUM_LIT:1> , self . n_t_controls ) <EOL> self . dYiqdCt = ( self . dPi_bar_rqdCt * self . sin_inner_t_theta ) . reshape ( - <NUM_LIT:1> , self . n_t_controls ) <EOL> self . dZiqdCt = ( self . dPi_bar_rqdCt * self . cos_inner_t_theta ) . reshape ( - <NUM_LIT:1> , self . n_t_controls ) <EOL> def copy ( self ) : <EOL> return copy . deepcopy ( self ) <EOL> def plot_geom ( self , ax , initial_color = '<STR_LIT:g>' , ffd_color = '<STR_LIT:k>' ) : <EOL> if initial_color : <EOL> ax . scatter ( self . Po [ : , <NUM_LIT:0> ] , self . Po [ : , <NUM_LIT:1> ] , c = initial_color , s = <NUM_LIT:50> , label = "<STR_LIT>" % self . name ) <EOL> ax . scatter ( self . Pi [ : , <NUM_LIT:0> ] , self . Pi [ : , <NUM_LIT:1> ] , c = initial_color , s = <NUM_LIT:50> ) <EOL> ax . plot ( self . Po [ : , <NUM_LIT:0> ] , self . Po [ : , <NUM_LIT:1> ] , c = initial_color ) <EOL> ax . plot ( self . Pi [ : , <NUM_LIT:0> ] , self . Pi [ : , <NUM_LIT:1> ] , c = initial_color ) <EOL> if ffd_color : <EOL> ax . scatter ( self . Po_bar [ : , <NUM_LIT:0> ] , self . Po_bar [ : , <NUM_LIT:1> ] , c = ffd_color , s = <NUM_LIT:50> , label = "<STR_LIT>" % self . name ) <EOL> ax . scatter ( self . Pi_bar [ : , <NUM_LIT:0> ] , self . Pi_bar [ : , <NUM_LIT:1> ] , c = ffd_color , s = <NUM_LIT:50> ) <EOL> ax . plot ( self . Po_bar [ : , <NUM_LIT:0> ] , self . Po_bar [ : , <NUM_LIT:1> ] , c = ffd_color ) <EOL> ax . plot ( self . Pi_bar [ : , <NUM_LIT:0> ] , self . Pi_bar [ : , <NUM_LIT:1> ] , c = ffd_color ) <EOL> def plot_centerline_spline ( self , ax , point_color = '<STR_LIT:r>' , line_color = '<STR_LIT:b>' ) : <EOL> ax . scatter ( self . Cc_bar [ : , <NUM_LIT:0> ] , self . Cc_bar [ : , <NUM_LIT:1> ] , c = point_color , s = <NUM_LIT:50> , label = "<STR_LIT>" % self . name ) <EOL> map_points = self . bsc_o ( np . linspace ( <NUM_LIT:0> , <NUM_LIT:1> , <NUM_LIT:100> ) ) <EOL> ax . plot ( map_points [ : , <NUM_LIT:0> ] , map_points [ : , <NUM_LIT:1> ] , label = "<STR_LIT>" , c = line_color ) <EOL> def plot_thickness_spline ( self , ax , point_color = '<STR_LIT:r>' , line_color = '<STR_LIT:b>' ) : <EOL> ax . scatter ( self . Ct_bar [ : , <NUM_LIT:0> ] , self . Ct_bar [ : , <NUM_LIT:1> ] , c = point_color , s = <NUM_LIT:50> , label = "<STR_LIT>" % self . name ) <EOL> map_points = self . bst_o ( np . linspace ( <NUM_LIT:0> , <NUM_LIT:1> , <NUM_LIT:100> ) ) <EOL> ax . plot ( map_points [ : , <NUM_LIT:0> ] , map_points [ : , <NUM_LIT:1> ] , label = "<STR_LIT>" , c = line_color ) <EOL> def deform ( self , delta_Cc , delta_Ct ) : <EOL> """<STR_LIT>""" <EOL> self . delta_Cc = delta_Cc <EOL> self . delta_Cc [ : , <NUM_LIT:0> ] *= self . x_mag <EOL> self . Cc_bar = self . Cc + self . delta_Cc <EOL> delta_Pc_o = self . bsc_o . calc ( self . Cc_bar ) <EOL> delta_Pc_i = self . bsc_i . calc ( self . Cc_bar ) <EOL> self . delta_Ct = delta_Ct <EOL> self . Ct_bar = self . Ct + self . delta_Ct <EOL> delta_Pt_o = self . bst_o . calc ( self . Ct_bar ) <EOL> delta_Pt_i = self . bst_i . calc ( self . Ct_bar ) <EOL> self . Po_bar = self . Po . copy ( ) <EOL> self . Pi_bar = self . Pi . copy ( ) <EOL> self . Po_bar [ : , <NUM_LIT:0> ] = delta_Pc_o [ : , <NUM_LIT:0> ] <EOL> self . Po_bar [ : , <NUM_LIT:1> ] = self . Po [ : , <NUM_LIT:1> ] + self . r_mag * ( delta_Pc_o [ : , <NUM_LIT:1> ] + delta_Pt_o [ : , <NUM_LIT:1> ] ) <EOL> self . Pi_bar [ : , <NUM_LIT:0> ] = delta_Pc_i [ : , <NUM_LIT:0> ] <EOL> self . Pi_bar [ : , <NUM_LIT:1> ] = self . Pi [ : , <NUM_LIT:1> ] + self . r_mag * ( delta_Pc_i [ : , <NUM_LIT:1> ] - delta_Pt_i [ : , <NUM_LIT:1> ] ) <EOL> self . outer_coords = Coordinates ( self . Po_bar , cartesian = False ) <EOL> self . inner_coords = Coordinates ( self . Pi_bar , cartesian = False ) <EOL> self . Po_bar_cart = self . outer_coords . cartesian <EOL> self . Xo = self . Po_bar_cart [ : , <NUM_LIT:0> ] <EOL> self . Yo = self . Po_bar_cart [ : , <NUM_LIT:1> ] <EOL> self . Zo = self . Po_bar_cart [ : , <NUM_LIT:2> ] <EOL> self . outer_stl . update_points ( self . Po_bar_cart ) <EOL> self . Pi_bar_cart = self . inner_coords . cartesian <EOL> self . Xi = self . Po_bar_cart [ : , <NUM_LIT:0> ] <EOL> self . Yi = self . Po_bar_cart [ : , <NUM_LIT:1> ] <EOL> self . Zi = self . Po_bar_cart [ : , <NUM_LIT:2> ] <EOL> self . inner_stl . update_points ( self . Pi_bar_cart ) <EOL> return self . Po_bar , self . Pi_bar <EOL> if __name__ == "<STR_LIT:__main__>" : <EOL> p = [ [ <NUM_LIT:0> , <NUM_LIT:0> , <NUM_LIT:0> ] , [ <NUM_LIT:0> , <NUM_LIT:0> , <NUM_LIT:1> ] , [ <NUM_LIT:0> , <NUM_LIT:1> , <NUM_LIT:0> ] ] <EOL> p_prime = Coordinates ( p , cartesian = True ) <EOL> print p_prime . cartesian <EOL> print p_prime . cylindrical <EOL> p = [ [ <NUM_LIT:0> , <NUM_LIT:0> , <NUM_LIT:0> ] , [ <NUM_LIT:0> , <NUM_LIT:1> , <NUM_LIT:0> ] , [ <NUM_LIT:0> , <NUM_LIT:1> , np . pi / <NUM_LIT:2> ] ] <EOL> p_prime = Coordinates ( p , cartesian = False ) <EOL> print p_prime . cartesian <EOL> print p_prime . cylindrical </s>
