Accueil | Cast3M

Télécharger J_el_TUB_CDAI_divers_chargements.dgibi
************************************************************************
************************************************************************
* Cas-test de validation pour le calcul de Jel en 2D et en 3D pour un
* Cylindre avec un Defaut Axisymetrique en peau Interne (CDAI) soumis a
* 3 types de chargements :
* - traction
* - pression sur les levres
* - temperature imposee
 
* Verification de la procedure G_THETA avec :
* - elements lineaires et quadratiques
* - le modele complet et le demi-modele
* - solution de type RESO et solution de type PASAPAS
* - un modele 2D et un modele 3D
 
* (Regroupement des cas-test :
* - rupt18.dgibi
* - rupt19.dgibi
* - rupt20.dgibi
* - rupt21.dgibi
* - rupt22.dgibi
* - rupt23.dgibi
* - rupt24.dgibi
* - rupt25.dgibi
* en un seul cas-test)
 
*   ^  z        ______
*   |          |      |
*   |          |      |
*   |          |      |
*   |          |      |
*   |          |      |
*   |          |      |
*   |          |      |
*   |          |      |
*   |---->     |--    |
*        r     |      |
*              |      |
*              |      |
*              |      |
*              |      |
*              |      |
*              |      |
*              |______|
 
OPTI 'ECHO' 0;
 
 
* -----------------------------------
*  DONNEES GEOMETRIQUES ET PHYSIQUES
* -----------------------------------
 
* H : HAUTEUR DU TUBE
* T : EPAISSEUR DU TUBE
* A : PROFONDEUR DE LA FISSURE
* RI : RAYON INTERNE
* RE : RAYON EXTERNE
H = 1.;
T = 60.E-3;  
A = T / 5; 
RI = T * 5;
RE = RI + T;
 
* PROPRIETE MATERIAUX
E0 = 185E3;
NU0 = 0.3;
ALFA0 = 13.08E-6;
 
* NIVEAU DE CHARGEMENT
P0T = -400.;
P0P = 400.;
DT0 = 300.;
 
 
* -----------------------
*  SOLUTIONS ANALYTIQUES
* -----------------------
 
* FONCTIONS D'INFLUENCE
I0 = 1.211;
I1 = 0.718;
 
* CONTRAINTES IMPOSEES POUR LE GRADIENT DE TEMPERATURE
SIG0 = ((E0*ALFA0*DT0)/(1-NU0)) * (RI/(3*T)) * ((2*(RE**2))/(RI*(RE+RI)) - 1);
SIG1 = -1. * ((E0*ALFA0*DT0)/(1-NU0));
 
* J ANALYTIQUES
JT = (1-(NU0**2)) * ((I0*(-1.*P0T)*((PI*A)**(1./2.)))**2) / E0;
JP = (1-(NU0**2)) * ((I0*(-1.*P0P)*((PI*A)**(1./2.)))**2) / E0;
JTH = (1-(NU0**2)) * ((((I0*SIG0)+ (I1*SIG1*(A/T)))*((PI*A)**(1./2.)))**2) / E0;
 
 
* -----------
*  MAILLAGE 
* -----------
 
* ON COMMENCE PAR UN MAILLAGE QUADRATIQUE 2D, ON LE TRANSFORMERA ENSUITE AU BESOIN
OPTI 'DIME' 2 'ELEM' 'QUA8' 'MODE' 'AXIS';
 
* 1ER CARRE AUTOUR DE LA POINTE DE FISSURE
NELP = 10;
DX1 = 200.E-6;
L1 = DX1 * NELP;
PF = (RI + A) 0.;
P1 = PF PLUS (L1 0.);
P2 = P1 PLUS (0. L1);
D12 = DROI NELP P1 P2;
CAR1 = D12 TRAN NELP ((0. - L1) 0.);
ELIM (CAR1 ET PF) 1.E-10;
CAR1 = COUL 'JAUN' CAR1;
 
* DEBUT DU 2EME CARRE AUTOUR DE LA POINTE DE FISSURE
DX2 = 400.E-6;
L2 = DX2 * NELP;
D34 = D12 HOMO (L2/L1) PF;
HOMO1 = REGL NELP D12 D34;
HOMO1 = COUL 'VERT' HOMO1;
CAR2 = CAR1 ET HOMO1;
 
* DERAFFINEMENT 3 POUR 1
LDERAF = DX2*4;
ELS = D34 ELEM (LECT 1 2 3);
P1 = ELS POIN 'INIT';
P2 = ELS POIN 'FINA';
P3 P4 = P2 P1 PLUS (LDERAF 0.);
COTE2 = DROI 1 P2 P3;
COTE3 = DROI 1 P3 P4;
COTE4 = DROI 1 P4 P1;
DERAF1 = DALL ELS COTE2 COTE3 COTE4;
 
* DERAFFINEMENT 3 POUR 1
ELS = D34 ELEM (LECT 4 5 6);
P1 = P2;
P2 = ELS POIN 'FINA';
P4 = P3;
P3 = P2 PLUS (LDERAF 0.);
COTE2 = DROI 1 P2 P3;
COTE3 = DROI 1 P3 P4;
COTE4 = DROI 1 P4 P1;
DERAF2 = DALL ELS COTE2 COTE3 COTE4;
 
* DERAFFINEMENT 4 POUR 2
ELS = D34 ELEM (LECT 7 8 9 10);
P1 = P2;
P2 = ELS POIN 'FINA';
P4 = P3;
P3 = P2 PLUS (LDERAF 0.);
COTE2 = DROI 1 P2 P3;
COTE3 = DROI 2 P3 P4;
COTE4 = DROI 1 P4 P1;
DERAF3 = DALL ELS COTE2 COTE3 COTE4;
 
* SYMETRIE PAR RAPPORT A LA DIAGONALE
DERAF0 = DERAF1 ET DERAF2 ET DERAF3;
SYME1 = (HOMO1 ET DERAF0) SYME 'DROI' PF (PF PLUS (1. 1.));
 
* COIN
COIN1 = COTE2 TRAN 1 (0. LDERAF);
MOITIE_DROITE = CAR2 ET DERAF0 ET SYME1 ET COIN1;
ELIM MOITIE_DROITE 1.E-10;
 
* SYMETRIE PAR RAPPORT A LA VERTICALE
MOITIE_GAUCHE = MOITIE_DROITE SYME 'DROI' PF (PF PLUS (0. 1.));
ELIM MOITIE_GAUCHE MOITIE_DROITE 1.E-10;
REC1 = MOITIE_DROITE ET MOITIE_GAUCHE;
 
* TRANSLATION JUSQU'A LA PEAU INTERNE
XMIN = MINI (COOR 1 REC1);
PGAUCHE = REC1 POIN 'DROI' (XMIN 0.) (XMIN 1.) 1.E-10;
CONREC1 = CONT REC1;
ELGAUCHE = CONREC1 ELEM 'APPUYE' PGAUCHE;
DX3 = 1.6E-3;
TRGAUCHE = ELGAUCHE TRAN ((RI - XMIN) 0.) 'DINI' DX3 'DFIN' (DX3*2);
 
* TRANSLATION JUSQU'A LA PEAU EXTERNE
XMAX = MAXI (COOR 1 REC1);
PDROITE = REC1 POIN 'DROI' (XMAX 0.) (XMAX 1.) 1.E-10;
ELDROITE = CONREC1 ELEM 'APPUYE' PDROITE;
TRDROITE = ELDROITE TRAN ((RE - XMAX) 0.) 'DINI' DX3 'DFIN' (DX3*2);
REC2 = REC1 ET TRDROITE ET TRGAUCHE;
 
* TRANSLATION JUSQU'AU FOND DU CYLINDRE
YMAX = MAXI (COOR 2 REC2);
PHAUT = REC2 POIN 'DROI' (0. YMAX) (1. YMAX) 1.E-10;
ELHAUT = (CONT REC2) ELEM 'APPUYE' PHAUT;
TRHAUT = ELHAUT TRAN (0. (((H/2) - YMAX))) 'DINI' DX3 'DFIN' T;
CYL_SUP = REC2 ET TRHAUT;
 
* SYMETRIE PAR RAPPORT A L'HORIZONTALE
CYL_INF = CYL_SUP SYME 'DROI' PF (PF PLUS (1. 0.));
 
* DEFINITION DU LIGAMENT
PEXT = CYL_SUP POIN 'PROC' (RE 0.);
CONSUP = CONT CYL_SUP;
PART1 = CONSUP ELEM 'COMPRIS' PF PEXT;
PART2 = CONSUP ELEM 'COMPRIS' PEXT PF;
LIGAMENT = PART1;
SI ((NBEL LIGAMENT) > (NBEL PART2));
    LIGAMENT = PART2;
FINSI;
ELIM LIGAMENT CYL_INF 1.E-10;
 
* DEFINITION DES LEVRES
PINT = CYL_SUP POIN 'PROC' (RI 0.);
PART1 = CONSUP ELEM 'COMPRIS' PF PINT;
PART2 = CONSUP ELEM 'COMPRIS' PINT PF;
LVSUP = PART1;
SI ((NBEL LVSUP) > (NBEL PART2));
    LVSUP = PART2;
FINSI;
SI (NEG (NOEU (LVSUP POIN 'INIT')) (NOEU PF));
    LVSUP = INVE LVSUP;
FINSI;
PINT = CYL_INF POIN 'PROC' (RI 0.);
CONINF = CONT CYL_INF;
PART1 = CONINF ELEM 'COMPRIS' PF PINT;
PART2 = CONINF ELEM 'COMPRIS' PINT PF;
LVINF = PART1;
SI ((NBEL LVINF) > (NBEL PART2));
    LVINF = PART2;
FINSI;
SI (NEG (NOEU (LVINF POIN 'FINA')) (NOEU PF));
    LVINF = INVE LVINF;
FINSI;
 
* BORDS DU HAUT ET DU BAS
YMAX = H/2;
YMIN = 0. - YMAX;
PHAUT = CYL_SUP POIN 'DROI' (0. YMAX) (1. YMAX) 1.E-10;
BHAUT = (CONT CYL_SUP) ELEM 'APPUYE' PHAUT;
SI ((COOR 1 (BHAUT POIN 'INIT')) < (COOR 1 (BHAUT POIN 'FINA')));
    BHAUT = INVE BHAUT;
FINSI;
PBAS = CYL_INF POIN 'DROI' (0. YMIN) (1. YMIN) 1.E-10;
BBAS = (CONT CYL_INF) ELEM 'APPUYE' PBAS;
SI ((COOR 1 (BBAS POIN 'INIT')) > (COOR 1 (BBAS POIN 'FINA')));
    BBAS = INVE BBAS;
FINSI;
 
 
* ---------------------
*  CALCUL PRELIMINAIRE 
* ---------------------
 
* ON FAIT UN PREMIER CALCUL AVEC RESO SUR LA STRUCTURE COMPLETE
* PUIS ON EN DEDUIT LES SOLUTIONS DANS LES AUTRES CAS SANS 
* FAIRE DE CALCULS SUPPLEMENTAIRES
 
* MODELE ET MATERIAU
OMEGA = CYL_SUP ET CYL_INF;
MO = MODE OMEGA 'MECANIQUE' 'ELASTIQUE';
MA = MATE MO 'YOUN' E0 'NU' NU0 'ALPH' ALFA0 'TALP' 0. 'TREF' 0.;
 
* CONDITIONS AUX LIMITES
CL1 = RELA 'ENSE' 'UZ' PHAUT;
CL2 = RELA 'ENSE' 'UZ' PBAS;
CL3 = BLOQ 'UZ' PF;
CL0 = CL1 ET CL2 ET CL3;
 
* CHARGEMENTS
MO_T = MODE (BHAUT ET BBAS) 'CHARGEMENT' 'PRESSION' 'CONS' 'DOMEGAF';
CH_T = PRES MO_T 'PRES' P0T;
F_T = BSIG MO_T CH_T;
 
LEVRES = LVSUP ET LVINF;
MO_P = MODE LEVRES 'CHARGEMENT' 'PRESSION' 'CONS' 'LEVRES';
CH_P = PRES MO_P 'PRES' P0P;
F_P = BSIG MO_P CH_P;
 
CHR = ((OMEGA COOR 1) - RI) / (RE - RI);
CHT0 = NOMC 'T' (DT0 * CHR);
SIG0 = THET MO MA CHT0;
F_TH = BSIG MO SIG0;
 
* APPEL A RESO
RI0 = RIGI MO MA;
UT UP UTH = RESO (RI0 ET CL0) F_T F_P F_TH;
U0 = MANU 'CHPO' OMEGA 2 'UR' 0. 'UZ' 0.;
 
 
* ------------
*  PROCEDURES 
* ------------
 
* PROCEDURE DE VERIFICATION DE LA PRECISION DU RESULTAT
DEBP CHECKERR ERR*'FLOTTANT' LIM*'FLOTTANT' M1*'MOT';
    MESS '  ERREUR EN' ' ' M1 ' = ' ERR '%';
    SI (ERR > LIM);
        MESS;
        ERRE 'ERREUR TROP GRANDE EN' ' ' M1;
    SINON;
        MESS '  => OK';
    FINSI;
FINP;
 
* PROCEDURE QUI APPELLE G_THETA
* DANS CETTE PROCEDURE, ON REDEFINIT LE MODELE, LE CHARGEMENT ET LES CONDITIONS
* AUX LIMITES. MEME SI ON N'APPELLE NI RESO, NI PASAPAS, ON EN A BESOIN POUR
* APPELER G_THETA
DEBP CASTEST DIMENSION*'MOT' INTERPOLATION*'MOT' SYSTEME*'MOT';
 
    DIM1 DIM2 = '2D' '3D';
    SI ((NEG DIMENSION DIM1) ET (NEG DIMENSION DIM2));
        ERR1 = CHAI 'DIMENSION DOIT VALOIR ''' DIM1 ''' ' 'OU' ' ''' DIM2 '';
        ERRE ERR1;
    FINSI;
    INT1 INT2 = 'LINEAIRE' 'QUADRATIQUE';
    SI ((NEG INTERPOLATION INT1) ET (NEG INTERPOLATION INT2));
        ERR1 = CHAI 'INTERPOLATION DOIT VALOIR ''' INT1 ''' ' 'OU' ' ''' INT2 '';
        ERRE ERR1;
    FINSI;
    MOD1 = MOT 'COMPLET';
    MOD2 = 'MOITIE';
    SI ((NEG SYSTEME MOD1) ET (NEG SYSTEME MOD2));
        ERR1 = CHAI 'SYSTEME DOIT VALOIR ''' MOD1 ''' ' 'OU' ' ''' MOD2 '';
        ERRE ERR1;
    FINSI;
 
*   DEFINITION DU SYSTEME
    SI (EGA SYSTEME 'COMPLET');
        OMEGA = CYL_SUP ET CYL_INF;
    SINON;
        OMEGA = CYL_SUP;
    FINSI;
 
*   PASSAGE EN LINEAIRE LE CAS ECHEANT
    SI (EGA INTERPOLATION 'LINEAIRE');
        OMEGA = CHAN 'LINE' OMEGA;
        LVSUP = CHAN 'LINE' LVSUP;
        LVINF = CHAN 'LINE' LVINF;
        BHAUT = CHAN 'LINE' BHAUT;
        SI (EGA SYSTEME 'MOITIE');
            LIGAMENT = CHAN 'LINE' LIGAMENT;
        SINON;
            BBAS = CHAN 'LINE' BBAS; 
        FINSI;
        SI (EGA DIMENSION '3D');
            PF = CHAN 'LINE' PF;
        FINSI;
    FINSI; 
 
*   MODELE ET MATERIAU
    MO = MODE OMEGA 'MECANIQUE' 'ELASTIQUE';
    MA = MATE MO 'YOUN' E0 'NU' NU0 'ALPH' ALFA0 'TALP' 0. 'TREF' 0.;
 
*   CONDITIONS AUX LIMITES
    SI (EGA DIMENSION '2D');
        U_AXE = 'UZ';
    SINON;
        U_AXE = 'UY';
    FINSI;
    PHAUT = CHAN 'POI1' BHAUT;
    CL1 = RELA 'ENSE' U_AXE PHAUT;
    CL0 = CL1 ET CL2;
    SI (EGA SYSTEME 'COMPLET');
        PBAS = CHAN 'POI1' BBAS;
        CL2A = RELA 'ENSE' U_AXE PBAS;
        CL2B = BLOQ U_AXE PF;
        CL2 = CL2A ET CL2B;
    SINON;
        CL2 = BLOQ U_AXE LIGAMENT;
    FINSI;
    CL0 = CL1 ET CL2;
    SI (EGA DIMENSION '3D'); 
        CL3 = OMEGA SYMT 'DEPL' PSYM11 PSYM12 PSYM13 1.E-10;
        CL4 = OMEGA SYMT 'DEPL' PSYM21 PSYM22 PSYM23 1.E-10;
        CL0 = CL0 ET CL3 ET CL4;
    FINSI;
 
*   CHARGEMENTS
    SI (EGA SYSTEME 'COMPLET');
        MO_T = MODE (BHAUT ET BBAS) 'CHARGEMENT' 'PRESSION' 'CONS' 'DOMEGAF';
    SINON;
        MO_T = MODE BHAUT 'CHARGEMENT' 'PRESSION' 'CONS' 'HAUT';
    FINSI;
    CH_T = PRES MO_T 'PRES' P0T;
    EVT = EVOL 'MANU' 'TEMP' (PROG 0. 1. 2. 3.) 'PRES' (PROG 0. 1. 0. 0.);
    CHA_T = CHAR 'PRES' CH_T EVT;
 
    LEVRES = LVSUP;
    SI (EGA SYSTEME 'COMPLET');
        LEVRES = LEVRES ET LVINF;
    FINSI;
    MO_P = MODE LEVRES 'CHARGEMENT' 'PRESSION' 'CONS' 'LEVRES';
    CH_P = PRES MO_P 'PRES' P0P;
    EVP = EVOL 'MANU' 'TEMP' (PROG 0. 1. 2. 3.) 'PRES' (PROG 0. 0. 1. 0.);
    CHA_P = CHAR 'PRES' CH_P EVP;
 
    CHR = ((OMEGA COOR 1) - RI) / (RE - RI);
    CHT0 = NOMC 'T' (DT0 * CHR);
    EVTH = EVOL 'MANU' (PROG 0. 1. 2. 3.) (PROG 0. 0. 0. 1.);
    CHA_TH = CHAR 'T' CHT0 EVTH;
 
*   SOLUTIONS A PARTIR DU CALCUL PRELIMINAIRE
    UT = REDU UT OMEGA;
    UP = REDU UP OMEGA;
    UTH = REDU UTH OMEGA;
 
*   APPELS A G_THETA
    J0 = TABL;
    J0.'OBJECTIF' = MOT 'J';
    J0.'FRONT_FISSURE' = PF;
    J0.'COUCHE' = 5;
    J0.'LEVRE_SUPERIEURE' = LVSUP;
    SI (EGA SYSTEME 'COMPLET');
        J0.'LEVRE_INFERIEURE' = LVINF;
    FINSI;
    J0.'CARACTERISTIQUES' = MA;
    J0.'BLOCAGES_MECANIQUES' = CL0;
*   TRACTION
    J = COPI J0;
    J.'MODELE' = MO ET MO_T;
    J.'PRESSION' = CH_T;
    J.'SOLUTION_RESO' = UT;
    G_THETA J;
    JELT = J.'RESULTATS';
*   PRESSION
    J = COPI J0;
    J.'MODELE' = MO ET MO_P;
    J.'PRESSION' = CH_P;
    J.'SOLUTION_RESO' = UP;
    G_THETA J;
    JELP = J.'RESULTATS';
*   TEMPERATURE
    J = COPI J0;
    J.'MODELE' = MO;
    J.'TEMPERATURES' = CHT0;
    J.'SOLUTION_RESO' = UTH;
    G_THETA J;
    JELTH = J.'RESULTATS';
 
*   ERREURS SUR J
    EMAXT EMAXP EMAXTH = 0.3 0.3 0.4;
    SI (EGA DIMENSION '3D');
        JELT = JELT.'GLOBAL';
        JELP = JELP.'GLOBAL';
        JELTH = JELTH.'GLOBAL';
    FINSI;        
    ERRT = ABS ((JELT-JT)/JT)*100.;
    ERRP = ABS ((JELP-JP)/JP)*100.;
    ERRTH = ABS ((JELTH-JTH)/JTH)*100.;
    MESS;
    MESS '- AVEC RESO :';
    CHECKERR ERRT EMAXT 'TRACTION';
    CHECKERR ERRP EMAXP 'PRESSION SUR LES LEVRES';
    CHECKERR ERRTH EMAXTH 'TEMPERATURE';
 
*   RESOLUTION AVEC PASAPAS DES 3 PROBLEMES (UN A CHAQUE PAS DE TEMPS)
*   AU PAS 1 : TRACTION SEULE
*   AU PAS 2 : PRESSION SUR LES LEVRES
*   AU PAS 3 : GRADIENT DE TEMPERATURE
*   ENCORE UNE FOIS, ON N'APPELLE PAS PASAPAS MAIS ON REMPLIT UNE TABLE ET
*   ON FAIT COMME SI ELLE ETAIT ISSUE DE PASAPAS
    CALCUL = TABL;
    CALCUL.'MODELE' = MO ET MO_T ET MO_P;
    CALCUL.'CARACTERISTIQUES' = MA;
    CALCUL.'BLOCAGES_MECANIQUES' = CL0;
    CALCUL.'CHARGEMENT' = CHA_T ET CHA_P ET CHA_TH;
    CALCUL.'TEMPS_CALCULES' = PROG 1. 2. 3.;
    TEM = TABL;
    TEM.(0) TEM.(1) TEM.(2) TEM.(3) = 0. 1. 2. 3.;
    CALCUL.'TEMPS' = TEM;
    USOL = TABL;
    USOL.(0) USOL.(1) USOL.(2) USOL.(3) = U0 UT UP UTH;
    CALCUL.'DEPLACEMENTS' = USOL;
    SSOL = TABL;
    REPE I 4;
        SSOL.(&I-1) = SIGM MO MA USOL.(&I-1);
    FIN I;
    SSOL.(3) = SSOL.(3) - (THET MO MA CHT0);
    CALCUL.'CONTRAINTES' = SSOL;
    CALCUL.'GRANDS_DEPLACEMENTS' = FAUX;
    CALCUL.'DYNAMIQUE' = FAUX;
 
*   APPEL A G_THETA
    J = COPI J0;
    OTER J 'BLOCAGES_MECANIQUES';
    OTER J 'CARACTERISTIQUES';
    J.'SOLUTION_PASAPAS' = CALCUL;
    G_THETA J;
    JELT = J.'RESULTATS'.(1);
    JELP = J.'RESULTATS'.(2);
    JELTH = J.'RESULTATS'.(3);
    SI (EGA DIMENSION '3D');
        JELT = JELT.'GLOBAL';
        JELP = JELP.'GLOBAL';
        JELTH = JELTH.'GLOBAL';
    FINSI;        
    ERRT = ABS ((JELT-JT)/JT)*100.;
    ERRP = ABS ((JELP-JP)/JP)*100.;
    ERRTH = ABS ((JELTH-JTH)/JTH)*100.;
 
*   ERREURS
    MESS;
    MESS '- AVEC PASAPAS :';
    CHECKERR ERRT EMAXT 'TRACTION';
    CHECKERR ERRP EMAXP 'PRESSION SUR LES LEVRES';
    CHECKERR ERRTH EMAXTH 'TEMPERATURE';
 
FINP;
 
 
* ------------
*  CALCULS 2D
* ------------
 
* ON LANCE LES TESTS
CASTEST '2D' 'QUADRATIQUE' 'COMPLET';
CASTEST '2D' 'QUADRATIQUE' 'MOITIE';
CASTEST '2D' 'LINEAIRE' 'COMPLET';
CASTEST '2D' 'LINEAIRE' 'MOITIE';
 
 
* ------------
*  CALCULS 3D
* ------------
 
* ASTUCE : PAS BESOIN DE RESOUDRE DE NOUVEAU LES 3 PROBLEMES MECANIQUES,
*          ON PEUT UTILISER LA SOLUTION 2D ET L'OPERATEUR 'MOIN' POUR
*          OBTENIR LA SOLUTION PAR DIFFERENCE ENTRE LES GEOMETRIES
*          INITIALES ET DEFORMEES
 
OMEGA_T = OMEGA PLUS UT;
OMEGA_P = OMEGA PLUS UP;
OMEGA_TH = OMEGA PLUS UTH;
 
OPTI 'DIME' 3 'ELEM' 'CU20';
 
* ANGLE DE ROTATION
ANG1 = 0.5;
* POINTS DE L'AXE DE ROTATION
PA = 0. 0. 0.;
PB = 0. 1. 0.;
* NOMBRE D'ELEMENTS DANS L'EPAISSEUR
NEP = 1;
 
OMEGA = OMEGA VOLU NEP 'ROTA' ANG1 PA PB;
OMEGA_T = OMEGA_T VOLU NEP 'ROTA' ANG1 PA PB;
OMEGA_P = OMEGA_P VOLU NEP 'ROTA' ANG1 PA PB;
OMEGA_TH = OMEGA_TH VOLU NEP 'ROTA' ANG1 PA PB;
 
* VOICI LES SOLUTIONS 3D :
UT = OMEGA_T MOIN OMEGA;
UP = OMEGA_P MOIN OMEGA;
UTH = OMEGA_TH MOIN OMEGA;
U0 = MANU 'CHPO' OMEGA 3 'UX' 0. 'UY' 0. 'UZ' 0.;
 
* IL FAUT REDEFINIR LES ENTITES GEOMETRIQUES IMPORTANTES
PSUP = (COOR 2 OMEGA) POIN 'EGSUPE' (0. - 1.E-10);
CYL_SUP = OMEGA ELEM 'APPUYE' PSUP;
PINF = (COOR 2 OMEGA) POIN 'EGINFE' 1.E-10;
CYL_INF = OMEGA ELEM 'APPUYE' PINF;
LVSUP = LVSUP ROTA NEP ANG1 PA PB;
LVINF = LVINF ROTA NEP ANG1 PA PB;
LIGAMENT = LIGAMENT ROTA NEP ANG1 PA PB;
BHAUT = BHAUT ROTA NEP ANG1 PA PB;
BBAS = BBAS ROTA NEP ANG1 PA PB;
PF = CERC NEP 'ROTA' ANG1 PF PA PB;
ELIM CYL_SUP (BHAUT ET LVSUP) 1.E-10;
ELIM CYL_INF (BBAS ET LVINF) 1.E-10;
ELIM OMEGA LIGAMENT 1.E-10;
ELIM PF LIGAMENT 1.E-10;
 
* ON AJOUTE LES POINTS POUR LES 2 PLANS DE SYMETRIE EN 3D
PSYM11 = CYL_SUP POIN 'PROC' (RI 0. 0.);
PSYM12 = CYL_SUP POIN 'PROC' (RE 0. 0.);
PSYM13 = CYL_SUP POIN 'PROC' (RI (H/2) 0.);
PSYM21 = CYL_SUP POIN 'PROC' (PSYM11 TOUR ANG1 PA PB);
PSYM22 = CYL_SUP POIN 'PROC' (PSYM12 TOUR ANG1 PA PB);
PSYM23 = CYL_SUP POIN 'PROC' (PSYM13 TOUR ANG1 PA PB);
 
* ON LANCE LES TESTS
CASTEST '3D' 'QUADRATIQUE' 'COMPLET';
CASTEST '3D' 'QUADRATIQUE' 'MOITIE';
CASTEST '3D' 'LINEAIRE' 'COMPLET';
CASTEST '3D' 'LINEAIRE' 'MOITIE';
 
FIN;