From: Pascale Noyret Date: Wed, 13 May 2009 08:36:55 +0000 (+0000) Subject: *** empty log message *** X-Git-Tag: V0_a1_SEPQT4~1 X-Git-Url: http://git.salome-platform.org/gitweb/?a=commitdiff_plain;h=32ceab9cca4d0d427b476c5aa5e9dffcfd4c11e1;p=modules%2Feficas.git *** empty log message *** --- diff --git a/Sous_epaisseur/sdistSEP.py b/Sous_epaisseur/sdistSEP.py index 9dec6f1b..abd21f2d 100644 --- a/Sous_epaisseur/sdistSEP.py +++ b/Sous_epaisseur/sdistSEP.py @@ -30,15 +30,13 @@ def main(): copyfiles('.',path_distrib,['LICENSE.TERMS','INSTALL','NEWS']) copyfiles('../Editeur',os.path.join(path_distrib,'Editeur'),['*.py','faqs.txt']) - copyfiles('../InterfaceQT4',os.path.join(path_distrib,'InterfaceQT'),['*.py','faqs.txt']) - copyfiles('../UiQT4',os.path.join(path_distrib,'Ui'),['*.ui','makefile']) - copyfiles('../SEP',os.path.join(path_distrib,'Openturns'),['*.py','*.ini']) + copyfiles('../InterfaceQT4',os.path.join(path_distrib,'InterfaceQT4'),['*.py','faqs.txt']) + copyfiles('../UiQT4',os.path.join(path_distrib,'UiQT4'),['*.ui','makefile']) + copyfiles('../Sous_epaisseur',os.path.join(path_distrib,'Sous_epaisseur'),['*.py','*.ini']) copyfiles('../Ihm',os.path.join(path_distrib,'Ihm'),['*.py']) copyfiles('../Extensions',os.path.join(path_distrib,'Extensions'),['*.py']) copyfiles('../Misc',os.path.join(path_distrib,'Misc'),['*.py']) copyfiles('../Accas',os.path.join(path_distrib,'Accas'),['*.py']) - copyfiles('../Accas',os.path.join(path_distrib,'Noyau'),['*.py']) - copyfiles('../Accas',os.path.join(path_distrib,'Validation'),['*.py']) # AIDE copyfiles('../AIDE',os.path.join(path_distrib,'AIDE'),['*.py']) copyfiles('../AIDE/fichiers',os.path.join(path_distrib,'AIDE','fichiers'),['*']) diff --git a/Sous_epaisseur/tube.py b/Sous_epaisseur/tube.py deleted file mode 100644 index 39f118d9..00000000 --- a/Sous_epaisseur/tube.py +++ /dev/null @@ -1,921 +0,0 @@ - - -# ========== -# PARAMETRES -# ========== - -# methode de filling: par les generatrices ou par des ecailles de tortue -methode = "tortue" # "generatrices" ou "tortue" - -# Longueur -longueur_tube = 6363.655 - -# Facteur d'amplification pour voir si le filling passe bien par les generatrices -amplification_factor = 1 - -# teste si le maillage gibi passe par tous les points -test_gibi = False -maillage_gibi = '/home/PROJETS/OUVERT/edf_anode/from_edf/2008_11_28/Tube/tube.mail.med' - -# MAILLAGE -# ======== - -type_maillage = "regle" # "regle" ou "libre" - -if type_maillage == "libre": - type_algo = "NETGEN_2D" # "BLSURF" ou "NETGEN_2D" - -if methode == "generatrices": - # nombre de segments dans la longueur des generatrices dans la zone de sous-epaisseur - nb_seg_generatrices = 150 # methode generatrices -else: - # nombre de segments dans la longueur des generatrices dans la zone de sous-epaisseur - nb_seg_generatrices = 5 # methode tortue - # distance entre 2 abscisses de points de mesure au dessous de laquelle on discretise avec nb_seg_petites_distances - # au lieu de nb_seg_generatrices - petite_distance = 100 - nb_seg_petites_distances = 3 - -# nombre de segments dans l'epaisseur du tube -nb_seg_ep= 3 - -# nombre de segments dans l'arc du tube entre 2 generatrices -if methode == "generatrices": - nb_seg_arc = 20 # methode generatrices: partition en 2 -else: - nb_seg_arc = 5 # methode tortue: partition en 8 - -# nombre de segments dans la longueur de transition -nb_seg_transition = 4 - -# nombre de segments dans la longueur d'amortissement -nb_seg_amortissement = 11 - -# =========== - -import salome -from geompy import * -import smesh - -import os, csv, math, pdb, time - -# DEBUT DES FONCTIONS -# =================================================== - -# supprime recursivement un element dans une liste -# @exemple : -# l=["a", " ", "b", " ", "c"] -# l = recursive_remove(l, " ") -# print l -# => ["a", "b", "c"] -def recursive_remove(l, txt): - finished = 0 - while not finished: - try: - l.remove(txt) - except ValueError: - finished = 1 - return l - - -## lit les valeurs des epaisseurs sur les generatrices a partir d'un fichier csv -# @return d_generatrices : dictionnaire dont la cle est le nom de la generatrice et la valeur est la liste des epaisseur sur cette generatrice -# @return l_abscisses: liste des abscisses des points de mesures -def read_generatrice(filename): - file = open(filename, "rb") - - reader = csv.reader(file, delimiter=';', quoting=csv.QUOTE_NONE) - - # Dictionnaire pour stocker les mesures sur une generatrice donnee - d_generatrices = {} - - # Liste des noms des generatrices - l_noms_generatrices = [] - - # Liste des abscisses des points de mesures - l_abscisses = [] - - for i, row in enumerate(reader): - # On supprime les cases vides - row = recursive_remove(row, "") - # Sur la premiere ligne, on releve les noms des generatrices - if i==0: - for nom in row: - if nom not in ["Abscisse", "Longueur"]: - # on initialise chaque generatrice avec une liste vide - d_generatrices[nom] = [] - l_noms_generatrices.append(nom) - # nombre de generatrices trouvees: - nb_generatrices = len(d_generatrices) - else: - # sur les lignes suivantes, on releve les mesures des epaisseurs - for j, nom_generatrice in enumerate(l_noms_generatrices): - # la liste des epaisseurs commence a 1 - epaisseur_str_fr = row[j+1] - # on convertit le format decimal francais par le format anglais - epaisseur_str = epaisseur_str_fr.replace(",", ".") - epaisseur = float(epaisseur_str) - d_generatrices[nom_generatrice].append(epaisseur) - # on ajoute la valeur de l'abscisse - abscisse_str = row[nb_generatrices+1] - abscisse = float(abscisse_str) - l_abscisses.append(abscisse) - - file.close() - - return d_generatrices, l_noms_generatrices, l_abscisses - -## lit les valeurs des angles definissant les generatrices a partir d'un fichier csv -# @return l_angles : liste des angles definissant les generatrices -def read_angles(filename): - file = open(filename, "rb") - - reader = csv.reader(file, delimiter=';', quoting=csv.QUOTE_NONE) - - # Liste des angles des generatrices - l_angles = [] - - for row in reader: - # On supprime les cases vides - row = recursive_remove(row, "") - - # si la ligne comporte 3 valeurs, on peut lire les angles - if len(row) == 3: - angle_str = row[2] - angle = float(angle_str) - l_angles.append(angle) - return l_angles - -## Cree une face a partir d'un nuage de points -# @param l_arcs_points : liste de points sur un contour -# @param closed_wire: flag pour savoir si le contour est deja ferme ou s'il faut le fermer -# @return face : une face passant par tous les points -# @warning: TODO: completer l'algo pour qu'il fonctionne si le nombre de generatrices est impair!!! -def MakeShellFromPoints(l_arcs_points, closed_wire = False): - # on cree les arcs pour chaque quart de cercle - nb_generatrices = len(l_arcs_points[0]) - nb_arcs_filling = nb_generatrices/2 - l_arcs_filling = [[] for i in range(nb_arcs_filling)] - if closed_wire: - if nb_generatrices%2 != 1: - raise "L'algo ne fonctionne pour l'instant qu'avec un nombre de generatrices impair" - else: - if nb_generatrices%2 != 0: - raise "L'algo ne fonctionne pour l'instant qu'avec un nombre de generatrices pair" - # Creation des arcs a partir des points - for arc_points in l_arcs_points: - if not closed_wire: - # Pour cloturer le contour - arc_points.append(arc_points[0]) - for i in range(nb_arcs_filling): - # 3 points a la meme abscisse sur 3 generatrices consecutives - arc = MakeArc(arc_points[2*i], arc_points[1+2*i], arc_points[2+2*i]) - l_arcs_filling[i].append(arc) - - - # on fait un filling pour tous les 1ers arcs, tous les 2emes arcs, ... jusqu'au 4e arc - l_faces_filling = [] - for i, arc_filling in enumerate(l_arcs_filling): - # On fait un filling entre 2 quarts de cercles - for j in range(0, len(arc_filling) - 1): - l_quart_arc = [arc_filling[j], arc_filling[j+1]] - compound_quart_arcs = MakeCompound(l_quart_arc) - quart_face_filling = MakeFilling(compound_quart_arcs, 0, 10, 1e-05, 1e-05, 0) - #addToStudy(quart_face_filling, "quart_face_filling_%i"%(j+1)) - l_faces_filling.append(quart_face_filling) - face = MakeShell(l_faces_filling) - return face - -## Cree une face a partir d'un nuage de points par un filling sur les generatrices -# @param l_arcs_points : liste de points sur un contour -# @param closed_wire: flag pour savoir si le contour est deja ferme ou s'il faut le fermer -# @return face : une face passant par tous les points -def MakeFillingFromPoints(l_arcs_points, closed_wire = False): - nb_generatrices = len(l_arcs_points[0]) - l_points_generatrices = [[] for i in range(nb_generatrices)] - l_generatrices = [] - # Creation des generatrices a partir des points - for arc_points in l_arcs_points: - for i, point in enumerate(arc_points): - # on ajoute le point dans la generatrice correspondante - l_points_generatrices[i].append(point) - for points_generatrice in l_points_generatrices: - generatrice_i = MakeInterpol(points_generatrice) - l_generatrices.append(generatrice_i) - if not closed_wire: - # Pour cloturer le contour - l_generatrices.append(l_generatrices[0]) - compound_generatrices = MakeCompound(l_generatrices) - face = MakeFilling(compound_generatrices, 0, 10, 1e-05, 1e-05, 0) - return face - - -# FIN DES FONCTIONS -# =============================================== - -time0 = time.time() - -# lecture des mesures dans le fichier csv -mesures_filename = os.path.join(dir_name, "mesure-transposee.csv") -d_generatrices, l_noms_generatrices, l_abscisses = read_generatrice(mesures_filename) - -# lecture des angles dans le fichier csv -angles_filename = os.path.join(dir_name, "mesure-angles.csv") -l_angles = read_angles(angles_filename) - -# dictionnaire indiquant les angles en fonction du nom de la generatrice -d_angles = {} -for nom, angle in zip(l_noms_generatrices, l_angles): - d_angles[nom] = angle - -time1 = time.time() - -print "Temps de lecture des fichiers: %.1f s."%(time1-time0) - -# Rq: pour conserver le point de plus faible epaisseur, il faut que la couture de la face de filling -# se situe sur la generatrice ou se situe ce point - -# Points et vecteurs de base - -P0 = MakeVertex(0, 0, 0) - -Vx = MakeVectorDXDYDZ(1000, 0, 0) -Vy = MakeVectorDXDYDZ(0, 1000, 0) -Vz = MakeVectorDXDYDZ(0, 0, 1) -plane_size = longueur_tube * 10 - -P_ext = MakeVertex(0, r_ext, 0) -cercle_ext = MakeRevolution(P_ext, Vx, 2*math.pi) -face_ext_milieu = MakePrismVecH(cercle_ext, Vx, l_abscisses[-1]) -addToStudy(face_ext_milieu, "face_ext_milieu") - -# initialisation de l'epaisseur minimale pour l'algo de recherche de l'epaisseur minimale -ep_min = r_ext -# initialisation de la liste de generatrices -l_generatrices = [] -# initialisation de la double liste d'arcs -l_arcs_points = [[] for abscisse in l_abscisses] -angle = 0 -# angle de la generatrice ou se situe l'epaisseur minimale -angle_ep_mini = 0 -# Point ou se situe l'epaisseur minimale -P_ep_mini = None -# indice de la generatrice ou se situe le point d'epaisseur minimale -i_generatrice_ep_mini = None -# Creation des generatrices -for i, nom_generatrice in enumerate(l_noms_generatrices): - angle += d_angles[nom_generatrice] - l_ep = d_generatrices[nom_generatrice] - l_points_ep = [] - j = 0 - for ep, abscisse in zip(l_ep, l_abscisses): - P_generatrice_tmp = MakeVertex(abscisse, r_ext - ep*amplification_factor, 0) - P_generatrice = MakeRotation(P_generatrice_tmp, Vx, math.radians(angle)) - #addToStudy(P_generatrice, "P_generatrice_%i"%(i+1)) - # pour la methode par les generatrices - l_points_ep.append(P_generatrice) - # pour la methode par les arcs - l_arcs_points[j].append(P_generatrice) - # Test sur l'epaisseur minimale - if ep < ep_min: - ep_min = ep - P_ep_mini = P_generatrice - i_generatrice_ep_mini = i - angle_ep_mini = angle - j+=1 - # creation des generatrices - generatrice = MakeInterpol(l_points_ep) - addToStudy(generatrice, "Generatrice_%s"%nom_generatrice) - l_generatrices.append(generatrice) - -print "epaisseur minimale mesuree: ", ep_min -addToStudy(P_ep_mini, "P_ep_mini") - -# On cree un objet contenant tous les points pour voir si la face generee passe par tous les points -l_points = [] -for arcs_points in l_arcs_points: - l_points += arcs_points -Tous_Points = MakeCompound(l_points) -addToStudy(Tous_Points, "Tous_Points") - -# methode par les generatrices -# ============================ - -# Pour s'assurer de passer par le point d'epaisseur minimale, -# on decalle la generatrice de recollement de la face et on l'ajoute a la fin pour fermer la face -l_generatrices = l_generatrices[i_generatrice_ep_mini:] + l_generatrices[:i_generatrice_ep_mini] + [l_generatrices[i_generatrice_ep_mini]] -# Creation de la face -compound_generatrices = MakeCompound(l_generatrices) -addToStudy(compound_generatrices, "compound_generatrices") -min_compound_generatrices = MinDistance(face_ext_milieu, compound_generatrices) -print "epaisseur minimale entre les generatrices et la face exterieure: ", min_compound_generatrices -if methode == "generatrices": - face_int_milieu_tmp = MakeFilling(compound_generatrices, 0, 10, 1e-05, 1e-05, 0) - face_int_milieu = ChangeOrientation(face_int_milieu_tmp) - addToStudy(face_int_milieu, "face_int_milieu") - min_distance = MinDistance(face_ext_milieu, face_int_milieu) - print "epaisseur minimale avec la methode generatrices: ", min_distance - - -# methode par les arcs avec filling arc par arc, 3 points par points (methode ecaille de tortue) -# ================================================================== - -if methode != "generatrices": - face_int_milieu = MakeShellFromPoints(l_arcs_points) - addToStudy(face_int_milieu, "face_int_milieu") - - min_distance = MinDistance(face_ext_milieu, face_int_milieu) - print "epaisseur minimale avec la methode ecaille de tortue: ", min_distance - -# => La face suit a la fois les generatrices et les arcs. L'epaisseur minimale est respectee - -# Partie complementaire de la face interieure du tube -# =================================================== - -# calcul de la longueur d'amortissement a partir de la formule de U4.PC.10-D -r_moyen = r_ext - ep_nominale/2. -l_amor_1 = 3/2.*math.sqrt(r_moyen**3/ep_nominale) -l_amor_2 = 3*math.sqrt(r_moyen*ep_nominale) -longueur_amortissement = max(l_amor_1, l_amor_2) - -print "longueur d'amortissement: ", longueur_amortissement - -# Longueur de transition entre tube deformé et longueur d'amortissement -longueur_transition = longueur_amortissement/5. -print "longueur de transition: ", longueur_transition - - -# On cree un nuage de points definissant les contours de la face de transition -r_int = r_ext - ep_nominale*amplification_factor -l_faces = [] -l_abscisses_transition = [] -# boucle pour traiter en meme temps le prolongement en bas et en haut -for i_abscisse, coef_translation in zip([0, -1], [-1, 1]): - l_arcs_points_transition = [] - - # On cree les points sur le cercle de la face de transition, en vis-a-vis des points du premier cercle de la face int - abscisse_transition = l_abscisses[i_abscisse] + coef_translation*longueur_transition - l_abscisses_transition.append(abscisse_transition) - l_arcs_points_transition_base = [] - angle = 0 - for nom_generatrice in l_noms_generatrices: - angle += d_angles[nom_generatrice] - P_transition_base_tmp = MakeVertex(abscisse_transition, r_int, 0) - P_transition_base = MakeRotation(P_transition_base_tmp, Vx, math.radians(angle)) - l_arcs_points_transition_base.append(P_transition_base) - - # contour bas - l_arcs_points_transition.append(l_arcs_points_transition_base) - # contour haut - l_arcs_points_transition.append(l_arcs_points[i_abscisse]) - if methode == "generatrices": - face_transition = MakeFillingFromPoints(l_arcs_points_transition) - if coef_translation == -1: - face_transition = ChangeOrientation(face_transition) - else: - face_transition = MakeShellFromPoints(l_arcs_points_transition) - addToStudy(face_transition, "face_transition") - l_faces.append(face_transition) - - - # On recupere le contour bas pour creer la face d'amortissement - - P_transition = MakeVertex(abscisse_transition, 0, 0) - l_edge_transition = GetShapesOnPlaneWithLocation(face_transition, ShapeType["EDGE"], Vx, P_transition, GEOM.ST_ON) - wire_bas_transition = MakeWire(l_edge_transition) - face_amortissement = MakePrismVecH(wire_bas_transition, Vx, coef_translation*longueur_amortissement) - addToStudy(face_amortissement, "face_amortissement") - l_faces.append(face_amortissement) - -l_faces.append(face_int_milieu) - -if methode == "generatrices": - face_int = MakeSewing(l_faces, 0.1) -else: - face_int = MakeShell(l_faces) -addToStudy(face_int, "face_int") - - -# Creation du tube solide -# ======================= - -# Face exterieure - -h_tube = l_abscisses[-1] - l_abscisses[0] + 2*(longueur_amortissement+longueur_transition) -abscisse_base_tube = l_abscisses[0]-(longueur_amortissement+longueur_transition) -cercle_ext_bas = MakeTranslation(cercle_ext, abscisse_base_tube, 0, 0) -face_ext_tmp = MakePrismVecH(cercle_ext_bas, Vx, h_tube) -# on tourne la face, pour ne pas avoir l'edge de couture -face_ext = MakeRotation(face_ext_tmp, Vx, math.radians(l_angles[0])) -#face_ext = face_ext_tmp -addToStudy(face_ext, "face_ext") - -# Face bas -P_bas_ext = MakeTranslation(P_ext, abscisse_base_tube, 0, 0) -cercle_int_bas = CreateGroup(face_int, ShapeType["EDGE"]) -l_cercle_int_bas = GetShapesOnPlaneWithLocation(face_int, ShapeType["EDGE"], Vx, P_bas_ext, GEOM.ST_ON) -UnionList(cercle_int_bas, l_cercle_int_bas) - -face_bas = MakeFaceWires([cercle_ext_bas, cercle_int_bas], 1) -addToStudy(face_bas, "face_bas") - -# Face haut -P_haut_ext = MakeTranslation(P_bas_ext, h_tube, 0, 0) -cercle_ext_haut = MakeTranslation(cercle_ext_bas, h_tube, 0, 0) -cercle_int_haut = CreateGroup(face_int, ShapeType["EDGE"]) -l_cercle_int_haut = GetShapesOnPlaneWithLocation(face_int, ShapeType["EDGE"], Vx, P_haut_ext, GEOM.ST_ON) -UnionList(cercle_int_haut, l_cercle_int_haut) - -face_haut = MakeFaceWires([cercle_ext_haut, cercle_int_haut], 1) -addToStudy(face_haut, "face_haut") - -l_faces_tube = [face_int, face_ext, face_bas, face_haut] -shell_tube = MakeShell(l_faces_tube) -addToStudy(shell_tube, "shell_tube") - -tube = MakeSolid([shell_tube]) -addToStudy(tube, "tube") - -time2 = time.time() - -print "Temps de creation de la geometrie: %.1f s."%(time2-time1) - - -# Partition pour que le maillage passe par les points de mesure -# ============================================================= - -l_plans_abscisses = [] - - -if methode == "generatrices": - l_abscisses_plan = [l_abscisses_transition[0]] + [l_abscisses[0], l_abscisses[-1]] + [l_abscisses_transition[1]] -else: - l_abscisses_plan = [l_abscisses_transition[0]] + l_abscisses[:] + [l_abscisses_transition[1]] - -# un plan par abscisse -for abscisse in l_abscisses_plan: - P_plan_part = MakeVertex(abscisse, 0, 0) - plan_part = MakePlane(P_plan_part, Vx, plane_size) - l_plans_abscisses.append(plan_part) - -P_axe_tube = MakeVertex(abscisse_base_tube, 0, 0) -axe_tube = MakePrismVecH(P_axe_tube, Vx, h_tube) - -# pour rabotter les plans de partition -P_bas = MakeVertex(abscisse_base_tube, 0, 0) -cylindre_int = MakeCylinder(P_bas, Vx, r_int/2., h_tube) - -angle = 0 -l_plans_diag = [] -# un plan sur toutes les generatrices -# Rq: si on cree un plan toutes les 2 generatrices uniquement, -# le maillage ne passera pas par les points de mesure de l'autre generatrice sur 2 -for i, nom_generatrice in enumerate(l_noms_generatrices): - angle += d_angles[nom_generatrice] - if (methode != "generatrices") or (methode == "generatrices" and (i==i_generatrice_ep_mini)): - # TODO: lorsque MakePartition fonctionnera mieux (!), - # supprimer le if ci-dessus pour toujours partitionner par les plans des generatrices - P_vec_plan_tmp = MakeVertex(0, r_int, 0) - P_vec_plan = MakeRotation(P_vec_plan_tmp, Vx, math.radians(angle)) - V_plan = MakeVector(P0, P_vec_plan) - plan_diag_tmp = MakePrismVecH(axe_tube, V_plan, 2*r_ext) - plan_diag = MakeCut(plan_diag_tmp, cylindre_int) - l_plans_diag.append(plan_diag) - -# TODO: lorsque MakePartition fonctionnera mieux (!), supprimer ce bloc -# car on aura partitionne par toutes les generatrices dans le bloc precedent. -if methode == "generatrices": - plan_oppose = MakeRotation(l_plans_diag[-1], Vx, math.pi) - l_plans_diag.append(plan_oppose) - -tous_plans = MakeCompound(l_plans_abscisses + l_plans_diag) -addToStudy(tous_plans, "tous_plans") - -plans_diag = MakeCompound(l_plans_diag) -addToStudy(plans_diag, "plans_diag") - -plans_abscisses = MakeCompound(l_plans_abscisses) -addToStudy(plans_abscisses, "plans_abscisses") - -if type_maillage == "regle": - # Partion d'un tube plein par la face_int - cylindre_tmp = MakeCylinder(P_bas, Vx, r_ext, h_tube) - # le cylindre ainsi cree a sa ligne de couture sur Z - # => on la decalle sur l'edge de couture de la face interieure: -pi/2 + angle_ep_mini - if methode == "generatrices": - cylindre = MakeRotation(cylindre_tmp, Vx, -math.pi/2. + math.radians(angle_ep_mini)) - else: - # en methode ecailles de tortue, la reparation plante si on tourne le cylindre - # mais reussi si le cylindre reste avec des edges de couture non paralleles! - cylindre = cylindre_tmp - addToStudy(cylindre, "cylindre") - - cylindre_part = MakePartition([cylindre, face_int]) - addToStudy(cylindre_part, "cylindre_part") - - # on recupere le solide correspondant au tube - P_tube = MakeVertex(abscisse_base_tube, (r_ext+r_int)/2., 0) - tube = GetBlockNearPoint(cylindre_part, P_tube) - addToStudy(tube, "tube") - - if methode == "generatrices": - # partition par plans diag puis plans abscisses - tube_part_tmp = MakePartition([tube], [plans_diag]) - addToStudy(tube_part_tmp, "tube_part_tmp") - tube_part = MakePartition([tube_part_tmp], [plans_abscisses]) - else: - # partition par plans abscisses puis plans diag - tube_part_tmp = MakePartition([tube], [plans_abscisses]) - addToStudy(tube_part_tmp, "tube_part_tmp") - tube_part = MakePartition([tube_part_tmp], [plans_diag]) - addToStudy(tube_part, "tube_part") - - tube_part_improved = CheckAndImprove(tube_part) - if not tube_part_improved: - print "pas de reparation effectuee" - tube_part_improved = tube_part - addToStudy(tube_part_improved, "tube_part_improved") - -else: - # on partitionne les faces du bas et du haut pour pouvoir imposer le nombre de segments dans l'epaisseur - face_bas_part = MakePartition([face_bas], [tous_plans], Limit=ShapeType["FACE"]) - face_haut_part = MakePartition([face_haut], [tous_plans], Limit=ShapeType["FACE"]) - # pour le maillage libre, on partitionne uniquement avec les points de mesure - # pour qu'ils soient contenus dans le maillage - if methode == "generatrices": - # pour la methode generatrices, il faut partitionner la face interieure - face_int_part = MakePartition([face_int], [tous_plans], Limit=ShapeType["FACE"]) - else: - # pour la methode tortue, il suffit de partitionner par les points - # (en fait, seuls manquent les points au milieu des arcs) - face_int_part = MakePartition([face_int], [Tous_Points], Limit=ShapeType["FACE"]) - - l_faces_tube = [face_int_part, face_ext, face_bas_part, face_haut_part] - shell_tube = MakeShell(l_faces_tube) - addToStudy(shell_tube, "shell_tube") - - tube_part_improved = MakeSolid([shell_tube]) - addToStudy(tube_part_improved, "tube_part_improved") - -time3 = time.time() - -print "Temps de partitionnement: %.1f s."%(time3-time2) - -# Sous-geometries pour les sous-maillages -# ======================================= - -# edges dans l'epaisseur -l_edges_bas = GetShapesOnPlaneWithLocation(tube_part_improved, ShapeType["EDGE"], Vx, P_bas, GEOM.ST_ON) -edges_bas = CreateGroup(tube_part_improved, ShapeType["EDGE"]) -UnionList(edges_bas, l_edges_bas) -edges_ep_bas = GetEdgesByLength (edges_bas, 0, r_ext-r_int + 1e-1) -addToStudyInFather(tube_part_improved, edges_ep_bas, "edges_ep_bas") - -if type_maillage == "libre": - # on recupere les faces bas et haut - l_faces_bas = GetShapesOnPlaneWithLocation(tube_part_improved, ShapeType["FACE"], Vx, P_bas, GEOM.ST_ON) - l_faces_haut = GetShapesOnPlaneWithLocation(tube_part_improved, ShapeType["FACE"], Vx, P_haut_ext, GEOM.ST_ON) - faces_extremites = CreateGroup(tube_part_improved, ShapeType["FACE"]) - UnionList(faces_extremites, l_faces_bas + l_faces_haut) - addToStudyInFather(tube_part_improved, faces_extremites, "faces_extremites") - -# edges sur les arcs -edges_arc_bas = CreateGroup(edges_bas, ShapeType["EDGE"]) -l_edges_ep_bas = SubShapeAllSorted(edges_ep_bas, ShapeType["EDGE"]) -UnionList(edges_arc_bas, l_edges_bas) -DifferenceList(edges_arc_bas, l_edges_ep_bas) -addToStudyInFather(tube_part_improved, edges_arc_bas, "edges_arc_bas") - -# on recupere la face interieure -l_face_int = GetShapesOnCylinder(tube_part_improved, ShapeType["FACE"], Vx, r_ext - ep_min, GEOM.ST_IN) -sub_face_int = CreateGroup(tube_part_improved, ShapeType["FACE"]) -UnionList(sub_face_int, l_face_int) -addToStudyInFather(tube_part_improved, sub_face_int, "SURF_INT") - -# on recupere les edges d'amortissement -P_bas_int = MakeTranslation(P_bas, 0, r_int, 0) -P_edge_tmp = MakeRotation(P_bas_int, Vx, math.radians(angle_ep_mini)) -P_edge_amortissement_1 = MakeTranslation(P_edge_tmp, longueur_amortissement/2., 0, 0) -P_edge_amortissement_2 = MakeTranslation(P_edge_tmp, h_tube-longueur_amortissement/2., 0, 0) -edge_amortissement_1 = GetEdgeNearPoint(tube_part_improved, P_edge_amortissement_1) -edge_amortissement_2 = GetEdgeNearPoint(tube_part_improved, P_edge_amortissement_2) -edges_amortissement = CreateGroup(tube_part_improved, ShapeType["EDGE"]) -UnionList(edges_amortissement, [edge_amortissement_1, edge_amortissement_2]) -addToStudyInFather(tube_part_improved, edges_amortissement, "edges_amortissement") - -# on recupere les edges de transition -P_edge_transition_1 = MakeTranslation(P_edge_tmp, longueur_transition/2. + longueur_amortissement, 0, 0) -P_edge_transition_2 = MakeTranslation(P_edge_tmp, h_tube-(longueur_transition/2. + longueur_amortissement), 0, 0) -edge_transition_1 = GetEdgeNearPoint(tube_part_improved, P_edge_transition_1) -edge_transition_2 = GetEdgeNearPoint(tube_part_improved, P_edge_transition_2) -edges_transition = CreateGroup(tube_part_improved, ShapeType["EDGE"]) -UnionList(edges_transition, [edge_transition_1, edge_transition_2]) -addToStudyInFather(tube_part_improved, edges_transition, "edges_transition") - -# on recupere les edges d'une generatrice - -axe_generatrice_tmp = MakeTranslation(Vx, 0, r_int, 0) -axe_generatrice = MakeRotation(axe_generatrice_tmp, Vx, math.radians(angle_ep_mini)) -l_edges_generatrice = GetShapesOnCylinder(tube_part_improved, ShapeType["EDGE"], axe_generatrice, ep_nominale, GEOM.ST_IN) -edges_generatrice = CreateGroup(tube_part_improved, ShapeType["EDGE"]) -UnionList(edges_generatrice, l_edges_generatrice) -DifferenceList(edges_generatrice, [edge_amortissement_1, edge_amortissement_2]+[edge_transition_1, edge_transition_2]) -addToStudyInFather(tube_part_improved, edges_generatrice, "edges_generatrice") - -# on recupere les edges d'une generatrice dont la longueur est inferieure a petite_distance -# pour pouvoir imposer moins de segments -if methode != "generatrices": - edges_petite_distance = CreateGroup(tube_part_improved, ShapeType["EDGE"]) - l_petite_distance = [] - for edge in l_edges_generatrice: - length = BasicProperties(edge)[0] - if length <= petite_distance: - l_petite_distance.append(edge) - UnionList(edges_petite_distance, l_petite_distance) - addToStudyInFather(tube_part_improved, edges_petite_distance, "edges_petite_distance") - - - -# Sous-geometries pour les groupes -# ================================ - -# on recupere la face interieure sans les zones de transition et d'amortissement -P_abs_first = MakeVertex(l_abscisses[0], 0, 0) -l_face_int_inf = GetShapesOnPlaneWithLocation(tube_part_improved, ShapeType["FACE"], Vx, P_abs_first, GEOM.ST_IN) -P_abs_last = MakeVertex(l_abscisses[-1], 0, 0) -l_face_int_sup = GetShapesOnPlaneWithLocation(tube_part_improved, ShapeType["FACE"], Vx, P_abs_last, GEOM.ST_OUT) -sub_face_int_milieu = CreateGroup(tube_part_improved, ShapeType["FACE"]) -UnionList(sub_face_int_milieu, l_face_int) -DifferenceList(sub_face_int_milieu, l_face_int_inf + l_face_int_sup) -addToStudyInFather(tube_part_improved, sub_face_int_milieu, "FaceIntM") - -# on recupere la face exterieure -l_face_ext = GetShapesOnCylinder(tube_part_improved, ShapeType["FACE"], Vx, r_ext, GEOM.ST_ON) -sub_face_ext = CreateGroup(tube_part_improved, ShapeType["FACE"]) -UnionList(sub_face_ext, l_face_ext) -addToStudyInFather(tube_part_improved, sub_face_ext, "SURF_EXT") - -# on recupere la face a l'extremite amont du tube -l_face_bas = GetShapesOnPlaneWithLocation(tube_part_improved, ShapeType["FACE"], Vx, P_bas, GEOM.ST_ON) -sub_face_bas = CreateGroup(tube_part_improved, ShapeType["FACE"]) -UnionList(sub_face_bas, l_face_bas) -addToStudyInFather(tube_part_improved, sub_face_bas, "CLGV") - -# on recupere la face a l'extremite aval du tube -l_face_haut = GetShapesOnPlaneWithLocation(tube_part_improved, ShapeType["FACE"], Vx, P_haut_ext, GEOM.ST_ON) -sub_face_haut = CreateGroup(tube_part_improved, ShapeType["FACE"]) -UnionList(sub_face_haut, l_face_haut) -addToStudyInFather(tube_part_improved, sub_face_haut, "EXTUBE") - -# On recupere les edges communs a face_int et face_haut -l_edge_int = GetShapesOnCylinderIDs(tube_part_improved, ShapeType["EDGE"], Vx, r_ext - ep_min, GEOM.ST_IN) -l_edge_haut = GetShapesOnPlaneWithLocationIDs(tube_part_improved, ShapeType["EDGE"], Vx, P_haut_ext, GEOM.ST_ON) -l_edge_bord_int_haut = [] -for id_edge in l_edge_int: - if id_edge in l_edge_haut: - l_edge_bord_int_haut.append(id_edge) - -edge_bord_int_haut = CreateGroup(tube_part_improved, ShapeType["EDGE"]) -UnionIDs(edge_bord_int_haut, l_edge_bord_int_haut) -addToStudyInFather(tube_part_improved, edge_bord_int_haut, "BORDTU") - -# on recupere le point d'epaisseur minimale -# avec la methode tortue le maillage passe forcement par le point d'epaisseur minimale -#if methode != "generatrices": -x, y, z = PointCoordinates(P_ep_mini) -P_ep_mini_sub = GetPoint(tube_part_improved, x, y, z, 1e-5) -addToStudyInFather(tube_part_improved, P_ep_mini_sub, "P_ep_mini") - -time4 = time.time() - -print "Temps de recuperation des sous-geometries: %.1f s."%(time4-time3) - -# MAILLAGE -# ======== - -# on divise le nombre de segments par le nombre de plans suivant les abscisses - -Maillage = smesh.Mesh(tube_part_improved, "Tube") - -if type_maillage == "regle": - algo1D = Maillage.Segment() - algo1D.NumberOfSegments(nb_seg_generatrices) - #algo1D.QuadraticMesh() - - Maillage.Quadrangle() - - Maillage.Hexahedron() -else: - # 30 - average_length = h_tube/(nb_seg_generatrices + 2*(nb_seg_amortissement+nb_seg_transition))/5. - # BLSURF est un algo 1D/2D - if type_algo == "BLSURF": - algo2D = Maillage.Triangle(algo=smesh.BLSURF) - algo2D.SetPhySize(average_length) - else: - algo1D = Maillage.Segment() - algo1D.LocalLength(average_length) - - Maillage.Triangle(smesh.NETGEN_2D) - - - #algo3D = Maillage.Tetrahedron(smesh.GHS3D) - algo3D = Maillage.Tetrahedron() - - -# hypotheses locales -# ================== - -# On impose finalement un maillage fin partout, seul moyen d'avoir plusieurs elements dans l'epaisseur -if type_maillage == "libre": - # 8 - average_length_extremites = average_length/8. - # BLSURF est un aglo 1D/2D - if type_algo == "BLSURF": - #algo2D = Maillage.Triangle(geom=faces_extremites, algo=smesh.BLSURF) - #algo2D.SetPhySize(average_length_extremites) - pass - else: - # 8 - algo1D = Maillage.Segment(faces_extremites) - algo1D.LocalLength(average_length_extremites) - -if type_maillage == "regle": - # dans l'epaisseur - algo1D = Maillage.Segment(edges_ep_bas) - algo1D.NumberOfSegments(nb_seg_ep) - algo1D.Propagation() - - # sur les arcs - algo1D = Maillage.Segment(edges_arc_bas) - algo1D.NumberOfSegments(nb_seg_arc) - algo1D.Propagation() - - # sur les longueurs d'amortissement - algo1D = Maillage.Segment(edges_amortissement) - algo1D.NumberOfSegments(nb_seg_amortissement) - algo1D.Propagation() - - # sur les longueurs de transition - algo1D = Maillage.Segment(edges_transition) - algo1D.NumberOfSegments(nb_seg_transition) - algo1D.Propagation() - - if methode == "tortue": - algo1D = Maillage.Segment(edges_petite_distance) - algo1D.NumberOfSegments(nb_seg_petites_distances) - algo1D.Propagation() - -Maillage.Compute() - -# on fait passer le maillage par le point de plus faible epaisseur -#if methode == "generatrices": - #x, y, z = PointCoordinates(P_ep_mini) - #id_node = Maillage.MeshToPassThroughAPoint(x, y, z) - ## on cree le groupe avec le point de plus faible epaisseur - #Maillage.MakeGroupByIds("P_ep_mini", smesh.NODE, [id_node]) - -# on a deja cree le groupe geometrique -#if methode != "generatrices": -Maillage.Group(P_ep_mini_sub) - -# conversion en quadratique (tres long à afficher) -#Maillage.ConvertToQuadratic(1) - -# on ajoute deux points aux extremites de l'axe du tube -x, y, z = PointCoordinates(P_bas) -id_p2 = Maillage.AddNode(x, y, z) -Maillage.MakeGroupByIds("P2", smesh.NODE, [id_p2]) - -id_p1 = Maillage.AddNode(x+h_tube, y, z) -Maillage.MakeGroupByIds("P1", smesh.NODE, [id_p1]) - -# Groupes -# ========= - -Maillage.Group(edge_bord_int_haut) -Maillage.Group(sub_face_int) -Maillage.Group(sub_face_int, "PEAUINT") -Maillage.Group(sub_face_ext) -Maillage.Group(sub_face_ext, "PEAUEXT") -Maillage.Group(sub_face_bas) -Maillage.Group(sub_face_bas, "FACE1") -Maillage.Group(sub_face_haut) -Maillage.Group(sub_face_haut, "FACE2") -group_nodes_int_milieu = Maillage.GroupOnGeom(sub_face_int_milieu, "NodesInt", smesh.NODE) -Maillage.GroupOnGeom(tube_part_improved, "VOL_TUBE", smesh.VOLUME) -Maillage.GroupOnGeom(tube_part_improved, "COUDE", smesh.VOLUME) - -if test_gibi: - time4 = time.time() - ([MaillageGibi], status) = smesh.smesh.CreateMeshesFromMED(maillage_gibi) - - # on met le maillage gibi dans le meme axe que le maillage SALOME - MaillageGibi.RotateObject(MaillageGibi, Vz, -math.pi/2., False) - - V_trans = MakeVectorDXDYDZ(-longueur_amortissement-longueur_transition+l_abscisses[0], 0, 0) - MaillageGibi.TranslateObject(MaillageGibi, V_trans, False) - - MaillageInt = MaillageGibi - - gibi_groupes = MaillageGibi.GetGroups() - - # on determine le groupe correspondant a la face interieure - group_int = None - for groupe in gibi_groupes: - name = groupe.GetName() - if name.strip() == "SURF_INT": - group_int = groupe - break - l_faces_int = group_int.GetIDs() - l_nodes_ids = [] - for face_id in l_faces_int: - l_nodes = MaillageGibi.GetElemNodes(face_id) - for node in l_nodes: - if node not in l_nodes_ids: - l_nodes_ids.append(node) - -time5 = time.time() - -print "Temps de generation du maillage: %.1f s."%(time5-time4) - -# Verifions si le maillage passe par les points de mesure -# ======================================================= - -# on cree un maillage avec les points de mesure -MaillageTousPoints = smesh.Mesh(Tous_Points, "Tous_Points") -## BUG: smesh ne peut pas creer un maillage de points! -#MaillageTousPoints.Compute() -# => On ajoute les points un par un... -l_tous_points = SubShapeAllSorted(Tous_Points, ShapeType["VERTEX"]) -for point in l_tous_points: - x, y, z = PointCoordinates(point) - MaillageTousPoints.AddNode(x, y, z) - -l_points_mesures_ids = MaillageTousPoints.GetNodesId() - -# on ajoute les noeuds mailles de la face interieure -if not test_gibi: - l_nodes_ids = group_nodes_int_milieu.GetIDs() - MaillageInt = Maillage -for node in l_nodes_ids: - # on recupere les coordonnees depuis le maillage global - x, y, z = MaillageInt.GetNodeXYZ(node) - # on ajoute ce noeud dans le maillage de points - MaillageTousPoints.AddNode(x, y, z) - -# on trouve les noeuds en double -tolerance = 1e0 -coincident_nodes = MaillageTousPoints.FindCoincidentNodes(tolerance) - -# nombre de points de mesure -nb_points = len(l_points) - -# nombre de noeuds en commun -nb_coincident_nodes = len(coincident_nodes) - -# nombre de points perdus -nb_points_perdus = nb_points - nb_coincident_nodes - -print "%i/%i points de mesure ont ete conserves dans le maillage"%(nb_coincident_nodes, nb_points) -if nb_points_perdus: - print "%i/%i points de mesure ont ete perdus dans le maillage"%(nb_points_perdus, nb_points) - -# affichage des points de mesure conserves -group_kept_nodes = MaillageTousPoints.CreateEmptyGroup(smesh.NODE, "Kept_measure_points") -l_id_coincident_nodes = [] -for l_id in coincident_nodes: - l_id_coincident_nodes += l_id -group_kept_nodes.Add(l_id_coincident_nodes) - -# affichage des points de mesure perdus -group_lost_nodes = MaillageTousPoints.CreateEmptyGroup(smesh.NODE, "Lost_measure_points") -l_id_lost_points = [] -for id_point in l_points_mesures_ids: - if id_point not in l_id_coincident_nodes: - l_id_lost_points.append(id_point) -group_lost_nodes.Add(l_id_lost_points) - -# On merge les noeuds en double -if coincident_nodes: - MaillageTousPoints.MergeNodes(coincident_nodes) - -# on regarde si le point d'epaisseur minimale fait partie des points garde -x_mini, y_mini, z_mini = PointCoordinates(P_ep_mini) -id_p_ep_mini = MaillageTousPoints.AddNode(x_mini, y_mini, z_mini) -MaillageTousPoints.MakeGroupByIds("P_ep_mini", smesh.NODE, [id_p_ep_mini]) - -coincident_nodes_ep_mini = MaillageTousPoints.FindCoincidentNodes(tolerance) -if coincident_nodes_ep_mini: - print "Le point d'epaisseur minimale a ete conserve" -else: - print "Le point d'epaisseur minimale a ete perdu" - -time6 = time.time() - -print "Temps de verification du maillage: %.1f s."%(time6-time5) - -print "Temps total: %.1f s."%(time6 - time0) - -salome.sg.updateObjBrowser(0) - diff --git a/convert/convert_SEP.py b/convert/convert_SEP.py new file mode 100644 index 00000000..b318cea5 --- /dev/null +++ b/convert/convert_SEP.py @@ -0,0 +1,39 @@ +# -*- coding: utf-8 -*- +# CONFIGURATION MANAGEMENT OF EDF VERSION +# ====================================================================== +# COPYRIGHT (C) 1991 - 2002 EDF R&D WWW.CODE-ASTER.ORG +# THIS PROGRAM IS FREE SOFTWARE; YOU CAN REDISTRIBUTE IT AND/OR MODIFY +# IT UNDER THE TERMS OF THE GNU GENERAL PUBLIC LICENSE AS PUBLISHED BY +# THE FREE SOFTWARE FOUNDATION; EITHER VERSION 2 OF THE LICENSE, OR +# (AT YOUR OPTION) ANY LATER VERSION. +# +# THIS PROGRAM IS DISTRIBUTED IN THE HOPE THAT IT WILL BE USEFUL, BUT +# WITHOUT ANY WARRANTY; WITHOUT EVEN THE IMPLIED WARRANTY OF +# MERCHANTABILITY OR FITNESS FOR A PARTICULAR PURPOSE. SEE THE GNU +# GENERAL PUBLIC LICENSE FOR MORE DETAILS. +# +# YOU SHOULD HAVE RECEIVED A COPY OF THE GNU GENERAL PUBLIC LICENSE +# ALONG WITH THIS PROGRAM; IF NOT, WRITE TO EDF R&D CODE_ASTER, +# 1 AVENUE DU GENERAL DE GAULLE, 92141 CLAMART CEDEX, FRANCE. +# +# +# ====================================================================== +""" +""" + +import parseur_python +from convert_python import * + +def entryPoint(): + """ + Retourne les informations nécessaires pour le chargeur de plugins + Ces informations sont retournées dans un dictionnaire + """ + return { + # Le nom du plugin + 'name' : 'SEP', + # La factory pour créer une instance du plugin + 'factory' : PythonParser, + } + + diff --git a/generator/generator_SEP.py b/generator/generator_SEP.py new file mode 100644 index 00000000..7439809a --- /dev/null +++ b/generator/generator_SEP.py @@ -0,0 +1,100 @@ +# -*- coding: utf-8 -*- +# CONFIGURATION MANAGEMENT OF EDF VERSION +# ====================================================================== +# COPYRIGHT (C) 1991 - 2002 EDF R&D WWW.CODE-ASTER.ORG +# THIS PROGRAM IS FREE SOFTWARE; YOU CAN REDISTRIBUTE IT AND/OR MODIFY +# IT UNDER THE TERMS OF THE GNU GENERAL PUBLIC LICENSE AS PUBLISHED BY +# THE FREE SOFTWARE FOUNDATION; EITHER VERSION 2 OF THE LICENSE, OR +# (AT YOUR OPTION) ANY LATER VERSION. +# +# THIS PROGRAM IS DISTRIBUTED IN THE HOPE THAT IT WILL BE USEFUL, BUT +# WITHOUT ANY WARRANTY; WITHOUT EVEN THE IMPLIED WARRANTY OF +# MERCHANTABILITY OR FITNESS FOR A PARTICULAR PURPOSE. SEE THE GNU +# GENERAL PUBLIC LICENSE FOR MORE DETAILS. +# +# YOU SHOULD HAVE RECEIVED A COPY OF THE GNU GENERAL PUBLIC LICENSE +# ALONG WITH THIS PROGRAM; IF NOT, WRITE TO EDF R&D CODE_ASTER, +# 1 AVENUE DU GENERAL DE GAULLE, 92141 CLAMART CEDEX, FRANCE. +# +# +# ====================================================================== +""" + Ce module contient le plugin generateur de fichier au format + SEP pour EFICAS. + +""" +import traceback +import types,string,re,os + +from generator_python import PythonGenerator + +def entryPoint(): + """ + Retourne les informations necessaires pour le chargeur de plugins + + Ces informations sont retournees dans un dictionnaire + """ + return { + # Le nom du plugin + 'name' : 'SEP', + # La factory pour creer une instance du plugin + 'factory' : SEPGenerator, + } + + +class SEPGenerator(PythonGenerator): + """ + Ce generateur parcourt un objet de type JDC et produit + un texte au format eficas et + un texte au format py + + """ + # Les extensions de fichier permis? + extensions=('.comm',) + + def gener(self,obj,format='brut'): + self.initDico() + self.text=PythonGenerator.gener(self,obj,format) + self.genereSEP() + return self.text + + def getTubePy(self) : + return self.texteTubePy + + def genereSEP(self) : + self.texteTubePy="# Parametres generes par Eficas \n" + for MC in self.dictMCVal.keys(): + ligne = MC +"="+ repr(self.dictMCVal[MC])+'\n' + self.texteTubePy=self.texteTubePy+ligne + print self.texteTubePy + + fichier=os.path.join(os.path.dirname(__file__),"tube.py") + f=open(fichier,'r') + for ligne in f.readlines(): + self.texteTubePy=self.texteTubePy+ligne + f.close + + def initDico(self) : + self.tube=0 + self.coude=0 + self.dictMCVal={} + self.texteTubePy="" + + def generMCSIMP(self,obj) : + """ + Convertit un objet MCSIMP en texte python + Remplit le dictionnaire des MCSIMP si nous ne sommes ni dans une loi, ni dans une variable + """ + s=PythonGenerator.generMCSIMP(self,obj) + self.dictMCVal[obj.nom]=obj.valeur + return s + + def generMACRO_ETAPE(self,obj): + print obj.nom + if obj.nom == "M_TUBE" : + self.tube=1 + if obj.nom == "M_COUDE" : + self.coude=1 + s=PythonGenerator.generMACRO_ETAPE(self,obj) + return s + diff --git a/generator/generator_openturns_study.py b/generator/generator_openturns_study.py index 8a76087b..1023eb90 100644 --- a/generator/generator_openturns_study.py +++ b/generator/generator_openturns_study.py @@ -87,20 +87,6 @@ class OpenturnsGenerator(PythonGenerator): self.dictTempo[obj.nom]=obj.valeur return s -# def generMCFACT(self,obj): -# # Il n est pas possible d utiliser obj.valeur qui n est pas -# # a jour pour les nouvelles variables ou les modifications -# if obj.nom in ( "Variables", "Files", ) : -# self.TraiteMCSIMP=0 -# self.dictTempo={} -# s=PythonGenerator.generMCFACT(self,obj) -# if obj.nom in ( "Variables", ) : -# self.listeVariables.append(self.dictTempo) -# self.dictTempo={} -# else : -# self.listeFichiers.append(self.dictTempo) -# self.TraiteMCSIMP=1 -# return s def generETAPE(self,obj): print "IDM: generETAPE dans generator_openturns_study.py" diff --git a/generator/tube.py b/generator/tube.py new file mode 100644 index 00000000..b353d8ff --- /dev/null +++ b/generator/tube.py @@ -0,0 +1,919 @@ + +# ========== +# PARAMETRES +# ========== + +# methode de filling: par les generatrices ou par des ecailles de tortue + +# Longueur +longueur_tube = 6363.655 + +# Facteur d'amplification pour voir si le filling passe bien par les generatrices +amplification_factor = 1 + +# teste si le maillage gibi passe par tous les points +test_gibi = False +maillage_gibi = '/home/PROJETS/OUVERT/edf_anode/from_edf/2008_11_28/Tube/tube.mail.med' + +# MAILLAGE +# ======== + +type_maillage = "regle" # "regle" ou "libre" + +if type_maillage == "libre": + type_algo = "NETGEN_2D" # "BLSURF" ou "NETGEN_2D" + +if methode == "generatrices": + # nombre de segments dans la longueur des generatrices dans la zone de sous-epaisseur + nb_seg_generatrices = 150 # methode generatrices +else: + # nombre de segments dans la longueur des generatrices dans la zone de sous-epaisseur + nb_seg_generatrices = 5 # methode tortue + # distance entre 2 abscisses de points de mesure au dessous de laquelle on discretise avec nb_seg_petites_distances + # au lieu de nb_seg_generatrices + petite_distance = 100 + nb_seg_petites_distances = 3 + +# nombre de segments dans l'epaisseur du tube +nb_seg_ep= 3 + +# nombre de segments dans l'arc du tube entre 2 generatrices +if methode == "generatrices": + nb_seg_arc = 20 # methode generatrices: partition en 2 +else: + nb_seg_arc = 5 # methode tortue: partition en 8 + +# nombre de segments dans la longueur de transition +nb_seg_transition = 4 + +# nombre de segments dans la longueur d'amortissement +nb_seg_amortissement = 11 + +# =========== + +import salome +from geompy import * +import smesh + +import os, csv, math, pdb, time + +# DEBUT DES FONCTIONS +# =================================================== + +# supprime recursivement un element dans une liste +# @exemple : +# l=["a", " ", "b", " ", "c"] +# l = recursive_remove(l, " ") +# print l +# => ["a", "b", "c"] +def recursive_remove(l, txt): + finished = 0 + while not finished: + try: + l.remove(txt) + except ValueError: + finished = 1 + return l + + +## lit les valeurs des epaisseurs sur les generatrices a partir d'un fichier csv +# @return d_generatrices : dictionnaire dont la cle est le nom de la generatrice et la valeur est la liste des epaisseur sur cette generatrice +# @return l_abscisses: liste des abscisses des points de mesures +def read_generatrice(filename): + file = open(filename, "rb") + + reader = csv.reader(file, delimiter=';', quoting=csv.QUOTE_NONE) + + # Dictionnaire pour stocker les mesures sur une generatrice donnee + d_generatrices = {} + + # Liste des noms des generatrices + l_noms_generatrices = [] + + # Liste des abscisses des points de mesures + l_abscisses = [] + + for i, row in enumerate(reader): + # On supprime les cases vides + row = recursive_remove(row, "") + # Sur la premiere ligne, on releve les noms des generatrices + if i==0: + for nom in row: + if nom not in ["Abscisse", "Longueur"]: + # on initialise chaque generatrice avec une liste vide + d_generatrices[nom] = [] + l_noms_generatrices.append(nom) + # nombre de generatrices trouvees: + nb_generatrices = len(d_generatrices) + else: + # sur les lignes suivantes, on releve les mesures des epaisseurs + for j, nom_generatrice in enumerate(l_noms_generatrices): + # la liste des epaisseurs commence a 1 + epaisseur_str_fr = row[j+1] + # on convertit le format decimal francais par le format anglais + epaisseur_str = epaisseur_str_fr.replace(",", ".") + epaisseur = float(epaisseur_str) + d_generatrices[nom_generatrice].append(epaisseur) + # on ajoute la valeur de l'abscisse + abscisse_str = row[nb_generatrices+1] + abscisse = float(abscisse_str) + l_abscisses.append(abscisse) + + file.close() + + return d_generatrices, l_noms_generatrices, l_abscisses + +## lit les valeurs des angles definissant les generatrices a partir d'un fichier csv +# @return l_angles : liste des angles definissant les generatrices +def read_angles(filename): + file = open(filename, "rb") + + reader = csv.reader(file, delimiter=';', quoting=csv.QUOTE_NONE) + + # Liste des angles des generatrices + l_angles = [] + + for row in reader: + # On supprime les cases vides + row = recursive_remove(row, "") + + # si la ligne comporte 3 valeurs, on peut lire les angles + if len(row) == 3: + angle_str = row[2] + angle = float(angle_str) + l_angles.append(angle) + return l_angles + +## Cree une face a partir d'un nuage de points +# @param l_arcs_points : liste de points sur un contour +# @param closed_wire: flag pour savoir si le contour est deja ferme ou s'il faut le fermer +# @return face : une face passant par tous les points +# @warning: TODO: completer l'algo pour qu'il fonctionne si le nombre de generatrices est impair!!! +def MakeShellFromPoints(l_arcs_points, closed_wire = False): + # on cree les arcs pour chaque quart de cercle + nb_generatrices = len(l_arcs_points[0]) + nb_arcs_filling = nb_generatrices/2 + l_arcs_filling = [[] for i in range(nb_arcs_filling)] + if closed_wire: + if nb_generatrices%2 != 1: + raise "L'algo ne fonctionne pour l'instant qu'avec un nombre de generatrices impair" + else: + if nb_generatrices%2 != 0: + raise "L'algo ne fonctionne pour l'instant qu'avec un nombre de generatrices pair" + # Creation des arcs a partir des points + for arc_points in l_arcs_points: + if not closed_wire: + # Pour cloturer le contour + arc_points.append(arc_points[0]) + for i in range(nb_arcs_filling): + # 3 points a la meme abscisse sur 3 generatrices consecutives + arc = MakeArc(arc_points[2*i], arc_points[1+2*i], arc_points[2+2*i]) + l_arcs_filling[i].append(arc) + + + # on fait un filling pour tous les 1ers arcs, tous les 2emes arcs, ... jusqu'au 4e arc + l_faces_filling = [] + for i, arc_filling in enumerate(l_arcs_filling): + # On fait un filling entre 2 quarts de cercles + for j in range(0, len(arc_filling) - 1): + l_quart_arc = [arc_filling[j], arc_filling[j+1]] + compound_quart_arcs = MakeCompound(l_quart_arc) + quart_face_filling = MakeFilling(compound_quart_arcs, 0, 10, 1e-05, 1e-05, 0) + #addToStudy(quart_face_filling, "quart_face_filling_%i"%(j+1)) + l_faces_filling.append(quart_face_filling) + face = MakeShell(l_faces_filling) + return face + +## Cree une face a partir d'un nuage de points par un filling sur les generatrices +# @param l_arcs_points : liste de points sur un contour +# @param closed_wire: flag pour savoir si le contour est deja ferme ou s'il faut le fermer +# @return face : une face passant par tous les points +def MakeFillingFromPoints(l_arcs_points, closed_wire = False): + nb_generatrices = len(l_arcs_points[0]) + l_points_generatrices = [[] for i in range(nb_generatrices)] + l_generatrices = [] + # Creation des generatrices a partir des points + for arc_points in l_arcs_points: + for i, point in enumerate(arc_points): + # on ajoute le point dans la generatrice correspondante + l_points_generatrices[i].append(point) + for points_generatrice in l_points_generatrices: + generatrice_i = MakeInterpol(points_generatrice) + l_generatrices.append(generatrice_i) + if not closed_wire: + # Pour cloturer le contour + l_generatrices.append(l_generatrices[0]) + compound_generatrices = MakeCompound(l_generatrices) + face = MakeFilling(compound_generatrices, 0, 10, 1e-05, 1e-05, 0) + return face + + +# FIN DES FONCTIONS +# =============================================== + +time0 = time.time() + +# lecture des mesures dans le fichier csv +mesures_filename = os.path.join(dir_name, "mesure-transposee.csv") +d_generatrices, l_noms_generatrices, l_abscisses = read_generatrice(mesures_filename) + +# lecture des angles dans le fichier csv +angles_filename = os.path.join(dir_name, "mesure-angles.csv") +l_angles = read_angles(angles_filename) + +# dictionnaire indiquant les angles en fonction du nom de la generatrice +d_angles = {} +for nom, angle in zip(l_noms_generatrices, l_angles): + d_angles[nom] = angle + +time1 = time.time() + +print "Temps de lecture des fichiers: %.1f s."%(time1-time0) + +# Rq: pour conserver le point de plus faible epaisseur, il faut que la couture de la face de filling +# se situe sur la generatrice ou se situe ce point + +# Points et vecteurs de base + +P0 = MakeVertex(0, 0, 0) + +Vx = MakeVectorDXDYDZ(1000, 0, 0) +Vy = MakeVectorDXDYDZ(0, 1000, 0) +Vz = MakeVectorDXDYDZ(0, 0, 1) +plane_size = longueur_tube * 10 + +P_ext = MakeVertex(0, r_ext, 0) +cercle_ext = MakeRevolution(P_ext, Vx, 2*math.pi) +face_ext_milieu = MakePrismVecH(cercle_ext, Vx, l_abscisses[-1]) +addToStudy(face_ext_milieu, "face_ext_milieu") + +# initialisation de l'epaisseur minimale pour l'algo de recherche de l'epaisseur minimale +ep_min = r_ext +# initialisation de la liste de generatrices +l_generatrices = [] +# initialisation de la double liste d'arcs +l_arcs_points = [[] for abscisse in l_abscisses] +angle = 0 +# angle de la generatrice ou se situe l'epaisseur minimale +angle_ep_mini = 0 +# Point ou se situe l'epaisseur minimale +P_ep_mini = None +# indice de la generatrice ou se situe le point d'epaisseur minimale +i_generatrice_ep_mini = None +# Creation des generatrices +for i, nom_generatrice in enumerate(l_noms_generatrices): + angle += d_angles[nom_generatrice] + l_ep = d_generatrices[nom_generatrice] + l_points_ep = [] + j = 0 + for ep, abscisse in zip(l_ep, l_abscisses): + P_generatrice_tmp = MakeVertex(abscisse, r_ext - ep*amplification_factor, 0) + P_generatrice = MakeRotation(P_generatrice_tmp, Vx, math.radians(angle)) + #addToStudy(P_generatrice, "P_generatrice_%i"%(i+1)) + # pour la methode par les generatrices + l_points_ep.append(P_generatrice) + # pour la methode par les arcs + l_arcs_points[j].append(P_generatrice) + # Test sur l'epaisseur minimale + if ep < ep_min: + ep_min = ep + P_ep_mini = P_generatrice + i_generatrice_ep_mini = i + angle_ep_mini = angle + j+=1 + # creation des generatrices + generatrice = MakeInterpol(l_points_ep) + addToStudy(generatrice, "Generatrice_%s"%nom_generatrice) + l_generatrices.append(generatrice) + +print "epaisseur minimale mesuree: ", ep_min +addToStudy(P_ep_mini, "P_ep_mini") + +# On cree un objet contenant tous les points pour voir si la face generee passe par tous les points +l_points = [] +for arcs_points in l_arcs_points: + l_points += arcs_points +Tous_Points = MakeCompound(l_points) +addToStudy(Tous_Points, "Tous_Points") + +# methode par les generatrices +# ============================ + +# Pour s'assurer de passer par le point d'epaisseur minimale, +# on decalle la generatrice de recollement de la face et on l'ajoute a la fin pour fermer la face +l_generatrices = l_generatrices[i_generatrice_ep_mini:] + l_generatrices[:i_generatrice_ep_mini] + [l_generatrices[i_generatrice_ep_mini]] +# Creation de la face +compound_generatrices = MakeCompound(l_generatrices) +addToStudy(compound_generatrices, "compound_generatrices") +min_compound_generatrices = MinDistance(face_ext_milieu, compound_generatrices) +print "epaisseur minimale entre les generatrices et la face exterieure: ", min_compound_generatrices +if methode == "generatrices": + face_int_milieu_tmp = MakeFilling(compound_generatrices, 0, 10, 1e-05, 1e-05, 0) + face_int_milieu = ChangeOrientation(face_int_milieu_tmp) + addToStudy(face_int_milieu, "face_int_milieu") + min_distance = MinDistance(face_ext_milieu, face_int_milieu) + print "epaisseur minimale avec la methode generatrices: ", min_distance + + +# methode par les arcs avec filling arc par arc, 3 points par points (methode ecaille de tortue) +# ================================================================== + +if methode != "generatrices": + face_int_milieu = MakeShellFromPoints(l_arcs_points) + addToStudy(face_int_milieu, "face_int_milieu") + + min_distance = MinDistance(face_ext_milieu, face_int_milieu) + print "epaisseur minimale avec la methode ecaille de tortue: ", min_distance + +# => La face suit a la fois les generatrices et les arcs. L'epaisseur minimale est respectee + +# Partie complementaire de la face interieure du tube +# =================================================== + +# calcul de la longueur d'amortissement a partir de la formule de U4.PC.10-D +r_moyen = r_ext - ep_nominale/2. +l_amor_1 = 3/2.*math.sqrt(r_moyen**3/ep_nominale) +l_amor_2 = 3*math.sqrt(r_moyen*ep_nominale) +longueur_amortissement = max(l_amor_1, l_amor_2) + +print "longueur d'amortissement: ", longueur_amortissement + +# Longueur de transition entre tube deformé et longueur d'amortissement +longueur_transition = longueur_amortissement/5. +print "longueur de transition: ", longueur_transition + + +# On cree un nuage de points definissant les contours de la face de transition +r_int = r_ext - ep_nominale*amplification_factor +l_faces = [] +l_abscisses_transition = [] +# boucle pour traiter en meme temps le prolongement en bas et en haut +for i_abscisse, coef_translation in zip([0, -1], [-1, 1]): + l_arcs_points_transition = [] + + # On cree les points sur le cercle de la face de transition, en vis-a-vis des points du premier cercle de la face int + abscisse_transition = l_abscisses[i_abscisse] + coef_translation*longueur_transition + l_abscisses_transition.append(abscisse_transition) + l_arcs_points_transition_base = [] + angle = 0 + for nom_generatrice in l_noms_generatrices: + angle += d_angles[nom_generatrice] + P_transition_base_tmp = MakeVertex(abscisse_transition, r_int, 0) + P_transition_base = MakeRotation(P_transition_base_tmp, Vx, math.radians(angle)) + l_arcs_points_transition_base.append(P_transition_base) + + # contour bas + l_arcs_points_transition.append(l_arcs_points_transition_base) + # contour haut + l_arcs_points_transition.append(l_arcs_points[i_abscisse]) + if methode == "generatrices": + face_transition = MakeFillingFromPoints(l_arcs_points_transition) + if coef_translation == -1: + face_transition = ChangeOrientation(face_transition) + else: + face_transition = MakeShellFromPoints(l_arcs_points_transition) + addToStudy(face_transition, "face_transition") + l_faces.append(face_transition) + + + # On recupere le contour bas pour creer la face d'amortissement + + P_transition = MakeVertex(abscisse_transition, 0, 0) + l_edge_transition = GetShapesOnPlaneWithLocation(face_transition, ShapeType["EDGE"], Vx, P_transition, GEOM.ST_ON) + wire_bas_transition = MakeWire(l_edge_transition) + face_amortissement = MakePrismVecH(wire_bas_transition, Vx, coef_translation*longueur_amortissement) + addToStudy(face_amortissement, "face_amortissement") + l_faces.append(face_amortissement) + +l_faces.append(face_int_milieu) + +if methode == "generatrices": + face_int = MakeSewing(l_faces, 0.1) +else: + face_int = MakeShell(l_faces) +addToStudy(face_int, "face_int") + + +# Creation du tube solide +# ======================= + +# Face exterieure + +h_tube = l_abscisses[-1] - l_abscisses[0] + 2*(longueur_amortissement+longueur_transition) +abscisse_base_tube = l_abscisses[0]-(longueur_amortissement+longueur_transition) +cercle_ext_bas = MakeTranslation(cercle_ext, abscisse_base_tube, 0, 0) +face_ext_tmp = MakePrismVecH(cercle_ext_bas, Vx, h_tube) +# on tourne la face, pour ne pas avoir l'edge de couture +face_ext = MakeRotation(face_ext_tmp, Vx, math.radians(l_angles[0])) +#face_ext = face_ext_tmp +addToStudy(face_ext, "face_ext") + +# Face bas +P_bas_ext = MakeTranslation(P_ext, abscisse_base_tube, 0, 0) +cercle_int_bas = CreateGroup(face_int, ShapeType["EDGE"]) +l_cercle_int_bas = GetShapesOnPlaneWithLocation(face_int, ShapeType["EDGE"], Vx, P_bas_ext, GEOM.ST_ON) +UnionList(cercle_int_bas, l_cercle_int_bas) + +face_bas = MakeFaceWires([cercle_ext_bas, cercle_int_bas], 1) +addToStudy(face_bas, "face_bas") + +# Face haut +P_haut_ext = MakeTranslation(P_bas_ext, h_tube, 0, 0) +cercle_ext_haut = MakeTranslation(cercle_ext_bas, h_tube, 0, 0) +cercle_int_haut = CreateGroup(face_int, ShapeType["EDGE"]) +l_cercle_int_haut = GetShapesOnPlaneWithLocation(face_int, ShapeType["EDGE"], Vx, P_haut_ext, GEOM.ST_ON) +UnionList(cercle_int_haut, l_cercle_int_haut) + +face_haut = MakeFaceWires([cercle_ext_haut, cercle_int_haut], 1) +addToStudy(face_haut, "face_haut") + +l_faces_tube = [face_int, face_ext, face_bas, face_haut] +shell_tube = MakeShell(l_faces_tube) +addToStudy(shell_tube, "shell_tube") + +tube = MakeSolid([shell_tube]) +addToStudy(tube, "tube") + +time2 = time.time() + +print "Temps de creation de la geometrie: %.1f s."%(time2-time1) + + +# Partition pour que le maillage passe par les points de mesure +# ============================================================= + +l_plans_abscisses = [] + + +if methode == "generatrices": + l_abscisses_plan = [l_abscisses_transition[0]] + [l_abscisses[0], l_abscisses[-1]] + [l_abscisses_transition[1]] +else: + l_abscisses_plan = [l_abscisses_transition[0]] + l_abscisses[:] + [l_abscisses_transition[1]] + +# un plan par abscisse +for abscisse in l_abscisses_plan: + P_plan_part = MakeVertex(abscisse, 0, 0) + plan_part = MakePlane(P_plan_part, Vx, plane_size) + l_plans_abscisses.append(plan_part) + +P_axe_tube = MakeVertex(abscisse_base_tube, 0, 0) +axe_tube = MakePrismVecH(P_axe_tube, Vx, h_tube) + +# pour rabotter les plans de partition +P_bas = MakeVertex(abscisse_base_tube, 0, 0) +cylindre_int = MakeCylinder(P_bas, Vx, r_int/2., h_tube) + +angle = 0 +l_plans_diag = [] +# un plan sur toutes les generatrices +# Rq: si on cree un plan toutes les 2 generatrices uniquement, +# le maillage ne passera pas par les points de mesure de l'autre generatrice sur 2 +for i, nom_generatrice in enumerate(l_noms_generatrices): + angle += d_angles[nom_generatrice] + if (methode != "generatrices") or (methode == "generatrices" and (i==i_generatrice_ep_mini)): + # TODO: lorsque MakePartition fonctionnera mieux (!), + # supprimer le if ci-dessus pour toujours partitionner par les plans des generatrices + P_vec_plan_tmp = MakeVertex(0, r_int, 0) + P_vec_plan = MakeRotation(P_vec_plan_tmp, Vx, math.radians(angle)) + V_plan = MakeVector(P0, P_vec_plan) + plan_diag_tmp = MakePrismVecH(axe_tube, V_plan, 2*r_ext) + plan_diag = MakeCut(plan_diag_tmp, cylindre_int) + l_plans_diag.append(plan_diag) + +# TODO: lorsque MakePartition fonctionnera mieux (!), supprimer ce bloc +# car on aura partitionne par toutes les generatrices dans le bloc precedent. +if methode == "generatrices": + plan_oppose = MakeRotation(l_plans_diag[-1], Vx, math.pi) + l_plans_diag.append(plan_oppose) + +tous_plans = MakeCompound(l_plans_abscisses + l_plans_diag) +addToStudy(tous_plans, "tous_plans") + +plans_diag = MakeCompound(l_plans_diag) +addToStudy(plans_diag, "plans_diag") + +plans_abscisses = MakeCompound(l_plans_abscisses) +addToStudy(plans_abscisses, "plans_abscisses") + +if type_maillage == "regle": + # Partion d'un tube plein par la face_int + cylindre_tmp = MakeCylinder(P_bas, Vx, r_ext, h_tube) + # le cylindre ainsi cree a sa ligne de couture sur Z + # => on la decalle sur l'edge de couture de la face interieure: -pi/2 + angle_ep_mini + if methode == "generatrices": + cylindre = MakeRotation(cylindre_tmp, Vx, -math.pi/2. + math.radians(angle_ep_mini)) + else: + # en methode ecailles de tortue, la reparation plante si on tourne le cylindre + # mais reussi si le cylindre reste avec des edges de couture non paralleles! + cylindre = cylindre_tmp + addToStudy(cylindre, "cylindre") + + cylindre_part = MakePartition([cylindre, face_int]) + addToStudy(cylindre_part, "cylindre_part") + + # on recupere le solide correspondant au tube + P_tube = MakeVertex(abscisse_base_tube, (r_ext+r_int)/2., 0) + tube = GetBlockNearPoint(cylindre_part, P_tube) + addToStudy(tube, "tube") + + if methode == "generatrices": + # partition par plans diag puis plans abscisses + tube_part_tmp = MakePartition([tube], [plans_diag]) + addToStudy(tube_part_tmp, "tube_part_tmp") + tube_part = MakePartition([tube_part_tmp], [plans_abscisses]) + else: + # partition par plans abscisses puis plans diag + tube_part_tmp = MakePartition([tube], [plans_abscisses]) + addToStudy(tube_part_tmp, "tube_part_tmp") + tube_part = MakePartition([tube_part_tmp], [plans_diag]) + addToStudy(tube_part, "tube_part") + + tube_part_improved = CheckAndImprove(tube_part) + if not tube_part_improved: + print "pas de reparation effectuee" + tube_part_improved = tube_part + addToStudy(tube_part_improved, "tube_part_improved") + +else: + # on partitionne les faces du bas et du haut pour pouvoir imposer le nombre de segments dans l'epaisseur + face_bas_part = MakePartition([face_bas], [tous_plans], Limit=ShapeType["FACE"]) + face_haut_part = MakePartition([face_haut], [tous_plans], Limit=ShapeType["FACE"]) + # pour le maillage libre, on partitionne uniquement avec les points de mesure + # pour qu'ils soient contenus dans le maillage + if methode == "generatrices": + # pour la methode generatrices, il faut partitionner la face interieure + face_int_part = MakePartition([face_int], [tous_plans], Limit=ShapeType["FACE"]) + else: + # pour la methode tortue, il suffit de partitionner par les points + # (en fait, seuls manquent les points au milieu des arcs) + face_int_part = MakePartition([face_int], [Tous_Points], Limit=ShapeType["FACE"]) + + l_faces_tube = [face_int_part, face_ext, face_bas_part, face_haut_part] + shell_tube = MakeShell(l_faces_tube) + addToStudy(shell_tube, "shell_tube") + + tube_part_improved = MakeSolid([shell_tube]) + addToStudy(tube_part_improved, "tube_part_improved") + +time3 = time.time() + +print "Temps de partitionnement: %.1f s."%(time3-time2) + +# Sous-geometries pour les sous-maillages +# ======================================= + +# edges dans l'epaisseur +l_edges_bas = GetShapesOnPlaneWithLocation(tube_part_improved, ShapeType["EDGE"], Vx, P_bas, GEOM.ST_ON) +edges_bas = CreateGroup(tube_part_improved, ShapeType["EDGE"]) +UnionList(edges_bas, l_edges_bas) +edges_ep_bas = GetEdgesByLength (edges_bas, 0, r_ext-r_int + 1e-1) +addToStudyInFather(tube_part_improved, edges_ep_bas, "edges_ep_bas") + +if type_maillage == "libre": + # on recupere les faces bas et haut + l_faces_bas = GetShapesOnPlaneWithLocation(tube_part_improved, ShapeType["FACE"], Vx, P_bas, GEOM.ST_ON) + l_faces_haut = GetShapesOnPlaneWithLocation(tube_part_improved, ShapeType["FACE"], Vx, P_haut_ext, GEOM.ST_ON) + faces_extremites = CreateGroup(tube_part_improved, ShapeType["FACE"]) + UnionList(faces_extremites, l_faces_bas + l_faces_haut) + addToStudyInFather(tube_part_improved, faces_extremites, "faces_extremites") + +# edges sur les arcs +edges_arc_bas = CreateGroup(edges_bas, ShapeType["EDGE"]) +l_edges_ep_bas = SubShapeAllSorted(edges_ep_bas, ShapeType["EDGE"]) +UnionList(edges_arc_bas, l_edges_bas) +DifferenceList(edges_arc_bas, l_edges_ep_bas) +addToStudyInFather(tube_part_improved, edges_arc_bas, "edges_arc_bas") + +# on recupere la face interieure +l_face_int = GetShapesOnCylinder(tube_part_improved, ShapeType["FACE"], Vx, r_ext - ep_min, GEOM.ST_IN) +sub_face_int = CreateGroup(tube_part_improved, ShapeType["FACE"]) +UnionList(sub_face_int, l_face_int) +addToStudyInFather(tube_part_improved, sub_face_int, "SURF_INT") + +# on recupere les edges d'amortissement +P_bas_int = MakeTranslation(P_bas, 0, r_int, 0) +P_edge_tmp = MakeRotation(P_bas_int, Vx, math.radians(angle_ep_mini)) +P_edge_amortissement_1 = MakeTranslation(P_edge_tmp, longueur_amortissement/2., 0, 0) +P_edge_amortissement_2 = MakeTranslation(P_edge_tmp, h_tube-longueur_amortissement/2., 0, 0) +edge_amortissement_1 = GetEdgeNearPoint(tube_part_improved, P_edge_amortissement_1) +edge_amortissement_2 = GetEdgeNearPoint(tube_part_improved, P_edge_amortissement_2) +edges_amortissement = CreateGroup(tube_part_improved, ShapeType["EDGE"]) +UnionList(edges_amortissement, [edge_amortissement_1, edge_amortissement_2]) +addToStudyInFather(tube_part_improved, edges_amortissement, "edges_amortissement") + +# on recupere les edges de transition +P_edge_transition_1 = MakeTranslation(P_edge_tmp, longueur_transition/2. + longueur_amortissement, 0, 0) +P_edge_transition_2 = MakeTranslation(P_edge_tmp, h_tube-(longueur_transition/2. + longueur_amortissement), 0, 0) +edge_transition_1 = GetEdgeNearPoint(tube_part_improved, P_edge_transition_1) +edge_transition_2 = GetEdgeNearPoint(tube_part_improved, P_edge_transition_2) +edges_transition = CreateGroup(tube_part_improved, ShapeType["EDGE"]) +UnionList(edges_transition, [edge_transition_1, edge_transition_2]) +addToStudyInFather(tube_part_improved, edges_transition, "edges_transition") + +# on recupere les edges d'une generatrice + +axe_generatrice_tmp = MakeTranslation(Vx, 0, r_int, 0) +axe_generatrice = MakeRotation(axe_generatrice_tmp, Vx, math.radians(angle_ep_mini)) +l_edges_generatrice = GetShapesOnCylinder(tube_part_improved, ShapeType["EDGE"], axe_generatrice, ep_nominale, GEOM.ST_IN) +edges_generatrice = CreateGroup(tube_part_improved, ShapeType["EDGE"]) +UnionList(edges_generatrice, l_edges_generatrice) +DifferenceList(edges_generatrice, [edge_amortissement_1, edge_amortissement_2]+[edge_transition_1, edge_transition_2]) +addToStudyInFather(tube_part_improved, edges_generatrice, "edges_generatrice") + +# on recupere les edges d'une generatrice dont la longueur est inferieure a petite_distance +# pour pouvoir imposer moins de segments +if methode != "generatrices": + edges_petite_distance = CreateGroup(tube_part_improved, ShapeType["EDGE"]) + l_petite_distance = [] + for edge in l_edges_generatrice: + length = BasicProperties(edge)[0] + if length <= petite_distance: + l_petite_distance.append(edge) + UnionList(edges_petite_distance, l_petite_distance) + addToStudyInFather(tube_part_improved, edges_petite_distance, "edges_petite_distance") + + + +# Sous-geometries pour les groupes +# ================================ + +# on recupere la face interieure sans les zones de transition et d'amortissement +P_abs_first = MakeVertex(l_abscisses[0], 0, 0) +l_face_int_inf = GetShapesOnPlaneWithLocation(tube_part_improved, ShapeType["FACE"], Vx, P_abs_first, GEOM.ST_IN) +P_abs_last = MakeVertex(l_abscisses[-1], 0, 0) +l_face_int_sup = GetShapesOnPlaneWithLocation(tube_part_improved, ShapeType["FACE"], Vx, P_abs_last, GEOM.ST_OUT) +sub_face_int_milieu = CreateGroup(tube_part_improved, ShapeType["FACE"]) +UnionList(sub_face_int_milieu, l_face_int) +DifferenceList(sub_face_int_milieu, l_face_int_inf + l_face_int_sup) +addToStudyInFather(tube_part_improved, sub_face_int_milieu, "FaceIntM") + +# on recupere la face exterieure +l_face_ext = GetShapesOnCylinder(tube_part_improved, ShapeType["FACE"], Vx, r_ext, GEOM.ST_ON) +sub_face_ext = CreateGroup(tube_part_improved, ShapeType["FACE"]) +UnionList(sub_face_ext, l_face_ext) +addToStudyInFather(tube_part_improved, sub_face_ext, "SURF_EXT") + +# on recupere la face a l'extremite amont du tube +l_face_bas = GetShapesOnPlaneWithLocation(tube_part_improved, ShapeType["FACE"], Vx, P_bas, GEOM.ST_ON) +sub_face_bas = CreateGroup(tube_part_improved, ShapeType["FACE"]) +UnionList(sub_face_bas, l_face_bas) +addToStudyInFather(tube_part_improved, sub_face_bas, "CLGV") + +# on recupere la face a l'extremite aval du tube +l_face_haut = GetShapesOnPlaneWithLocation(tube_part_improved, ShapeType["FACE"], Vx, P_haut_ext, GEOM.ST_ON) +sub_face_haut = CreateGroup(tube_part_improved, ShapeType["FACE"]) +UnionList(sub_face_haut, l_face_haut) +addToStudyInFather(tube_part_improved, sub_face_haut, "EXTUBE") + +# On recupere les edges communs a face_int et face_haut +l_edge_int = GetShapesOnCylinderIDs(tube_part_improved, ShapeType["EDGE"], Vx, r_ext - ep_min, GEOM.ST_IN) +l_edge_haut = GetShapesOnPlaneWithLocationIDs(tube_part_improved, ShapeType["EDGE"], Vx, P_haut_ext, GEOM.ST_ON) +l_edge_bord_int_haut = [] +for id_edge in l_edge_int: + if id_edge in l_edge_haut: + l_edge_bord_int_haut.append(id_edge) + +edge_bord_int_haut = CreateGroup(tube_part_improved, ShapeType["EDGE"]) +UnionIDs(edge_bord_int_haut, l_edge_bord_int_haut) +addToStudyInFather(tube_part_improved, edge_bord_int_haut, "BORDTU") + +# on recupere le point d'epaisseur minimale +# avec la methode tortue le maillage passe forcement par le point d'epaisseur minimale +#if methode != "generatrices": +x, y, z = PointCoordinates(P_ep_mini) +P_ep_mini_sub = GetPoint(tube_part_improved, x, y, z, 1e-5) +addToStudyInFather(tube_part_improved, P_ep_mini_sub, "P_ep_mini") + +time4 = time.time() + +print "Temps de recuperation des sous-geometries: %.1f s."%(time4-time3) + +# MAILLAGE +# ======== + +# on divise le nombre de segments par le nombre de plans suivant les abscisses + +Maillage = smesh.Mesh(tube_part_improved, "Tube") + +if type_maillage == "regle": + algo1D = Maillage.Segment() + algo1D.NumberOfSegments(nb_seg_generatrices) + #algo1D.QuadraticMesh() + + Maillage.Quadrangle() + + Maillage.Hexahedron() +else: + # 30 + average_length = h_tube/(nb_seg_generatrices + 2*(nb_seg_amortissement+nb_seg_transition))/5. + # BLSURF est un algo 1D/2D + if type_algo == "BLSURF": + algo2D = Maillage.Triangle(algo=smesh.BLSURF) + algo2D.SetPhySize(average_length) + else: + algo1D = Maillage.Segment() + algo1D.LocalLength(average_length) + + Maillage.Triangle(smesh.NETGEN_2D) + + + #algo3D = Maillage.Tetrahedron(smesh.GHS3D) + algo3D = Maillage.Tetrahedron() + + +# hypotheses locales +# ================== + +# On impose finalement un maillage fin partout, seul moyen d'avoir plusieurs elements dans l'epaisseur +if type_maillage == "libre": + # 8 + average_length_extremites = average_length/8. + # BLSURF est un aglo 1D/2D + if type_algo == "BLSURF": + #algo2D = Maillage.Triangle(geom=faces_extremites, algo=smesh.BLSURF) + #algo2D.SetPhySize(average_length_extremites) + pass + else: + # 8 + algo1D = Maillage.Segment(faces_extremites) + algo1D.LocalLength(average_length_extremites) + +if type_maillage == "regle": + # dans l'epaisseur + algo1D = Maillage.Segment(edges_ep_bas) + algo1D.NumberOfSegments(nb_seg_ep) + algo1D.Propagation() + + # sur les arcs + algo1D = Maillage.Segment(edges_arc_bas) + algo1D.NumberOfSegments(nb_seg_arc) + algo1D.Propagation() + + # sur les longueurs d'amortissement + algo1D = Maillage.Segment(edges_amortissement) + algo1D.NumberOfSegments(nb_seg_amortissement) + algo1D.Propagation() + + # sur les longueurs de transition + algo1D = Maillage.Segment(edges_transition) + algo1D.NumberOfSegments(nb_seg_transition) + algo1D.Propagation() + + if methode == "tortue": + algo1D = Maillage.Segment(edges_petite_distance) + algo1D.NumberOfSegments(nb_seg_petites_distances) + algo1D.Propagation() + +Maillage.Compute() + +# on fait passer le maillage par le point de plus faible epaisseur +#if methode == "generatrices": + #x, y, z = PointCoordinates(P_ep_mini) + #id_node = Maillage.MeshToPassThroughAPoint(x, y, z) + ## on cree le groupe avec le point de plus faible epaisseur + #Maillage.MakeGroupByIds("P_ep_mini", smesh.NODE, [id_node]) + +# on a deja cree le groupe geometrique +#if methode != "generatrices": +Maillage.Group(P_ep_mini_sub) + +# conversion en quadratique (tres long à afficher) +#Maillage.ConvertToQuadratic(1) + +# on ajoute deux points aux extremites de l'axe du tube +x, y, z = PointCoordinates(P_bas) +id_p2 = Maillage.AddNode(x, y, z) +Maillage.MakeGroupByIds("P2", smesh.NODE, [id_p2]) + +id_p1 = Maillage.AddNode(x+h_tube, y, z) +Maillage.MakeGroupByIds("P1", smesh.NODE, [id_p1]) + +# Groupes +# ========= + +Maillage.Group(edge_bord_int_haut) +Maillage.Group(sub_face_int) +Maillage.Group(sub_face_int, "PEAUINT") +Maillage.Group(sub_face_ext) +Maillage.Group(sub_face_ext, "PEAUEXT") +Maillage.Group(sub_face_bas) +Maillage.Group(sub_face_bas, "FACE1") +Maillage.Group(sub_face_haut) +Maillage.Group(sub_face_haut, "FACE2") +group_nodes_int_milieu = Maillage.GroupOnGeom(sub_face_int_milieu, "NodesInt", smesh.NODE) +Maillage.GroupOnGeom(tube_part_improved, "VOL_TUBE", smesh.VOLUME) +Maillage.GroupOnGeom(tube_part_improved, "COUDE", smesh.VOLUME) + +if test_gibi: + time4 = time.time() + ([MaillageGibi], status) = smesh.smesh.CreateMeshesFromMED(maillage_gibi) + + # on met le maillage gibi dans le meme axe que le maillage SALOME + MaillageGibi.RotateObject(MaillageGibi, Vz, -math.pi/2., False) + + V_trans = MakeVectorDXDYDZ(-longueur_amortissement-longueur_transition+l_abscisses[0], 0, 0) + MaillageGibi.TranslateObject(MaillageGibi, V_trans, False) + + MaillageInt = MaillageGibi + + gibi_groupes = MaillageGibi.GetGroups() + + # on determine le groupe correspondant a la face interieure + group_int = None + for groupe in gibi_groupes: + name = groupe.GetName() + if name.strip() == "SURF_INT": + group_int = groupe + break + l_faces_int = group_int.GetIDs() + l_nodes_ids = [] + for face_id in l_faces_int: + l_nodes = MaillageGibi.GetElemNodes(face_id) + for node in l_nodes: + if node not in l_nodes_ids: + l_nodes_ids.append(node) + +time5 = time.time() + +print "Temps de generation du maillage: %.1f s."%(time5-time4) + +# Verifions si le maillage passe par les points de mesure +# ======================================================= + +# on cree un maillage avec les points de mesure +MaillageTousPoints = smesh.Mesh(Tous_Points, "Tous_Points") +## BUG: smesh ne peut pas creer un maillage de points! +#MaillageTousPoints.Compute() +# => On ajoute les points un par un... +l_tous_points = SubShapeAllSorted(Tous_Points, ShapeType["VERTEX"]) +for point in l_tous_points: + x, y, z = PointCoordinates(point) + MaillageTousPoints.AddNode(x, y, z) + +l_points_mesures_ids = MaillageTousPoints.GetNodesId() + +# on ajoute les noeuds mailles de la face interieure +if not test_gibi: + l_nodes_ids = group_nodes_int_milieu.GetIDs() + MaillageInt = Maillage +for node in l_nodes_ids: + # on recupere les coordonnees depuis le maillage global + x, y, z = MaillageInt.GetNodeXYZ(node) + # on ajoute ce noeud dans le maillage de points + MaillageTousPoints.AddNode(x, y, z) + +# on trouve les noeuds en double +tolerance = 1e0 +coincident_nodes = MaillageTousPoints.FindCoincidentNodes(tolerance) + +# nombre de points de mesure +nb_points = len(l_points) + +# nombre de noeuds en commun +nb_coincident_nodes = len(coincident_nodes) + +# nombre de points perdus +nb_points_perdus = nb_points - nb_coincident_nodes + +print "%i/%i points de mesure ont ete conserves dans le maillage"%(nb_coincident_nodes, nb_points) +if nb_points_perdus: + print "%i/%i points de mesure ont ete perdus dans le maillage"%(nb_points_perdus, nb_points) + +# affichage des points de mesure conserves +group_kept_nodes = MaillageTousPoints.CreateEmptyGroup(smesh.NODE, "Kept_measure_points") +l_id_coincident_nodes = [] +for l_id in coincident_nodes: + l_id_coincident_nodes += l_id +group_kept_nodes.Add(l_id_coincident_nodes) + +# affichage des points de mesure perdus +group_lost_nodes = MaillageTousPoints.CreateEmptyGroup(smesh.NODE, "Lost_measure_points") +l_id_lost_points = [] +for id_point in l_points_mesures_ids: + if id_point not in l_id_coincident_nodes: + l_id_lost_points.append(id_point) +group_lost_nodes.Add(l_id_lost_points) + +# On merge les noeuds en double +if coincident_nodes: + MaillageTousPoints.MergeNodes(coincident_nodes) + +# on regarde si le point d'epaisseur minimale fait partie des points garde +x_mini, y_mini, z_mini = PointCoordinates(P_ep_mini) +id_p_ep_mini = MaillageTousPoints.AddNode(x_mini, y_mini, z_mini) +MaillageTousPoints.MakeGroupByIds("P_ep_mini", smesh.NODE, [id_p_ep_mini]) + +coincident_nodes_ep_mini = MaillageTousPoints.FindCoincidentNodes(tolerance) +if coincident_nodes_ep_mini: + print "Le point d'epaisseur minimale a ete conserve" +else: + print "Le point d'epaisseur minimale a ete perdu" + +time6 = time.time() + +print "Temps de verification du maillage: %.1f s."%(time6-time5) + +print "Temps total: %.1f s."%(time6 - time0) + +salome.sg.updateObjBrowser(0) +