* fichier : g_theta_utilisateur_2.dgibi
****************************************************
* Section : Mecanique rupture
****************************************************
*                                                  *
*      VERIFICATION DE LA PROCEDURE G_THETA        *
*    POUR LE CALCUL DE G POUR UNE FISSURE PLANE    *
*           A FOND DROIT DANS UN CUBE              *
*                        -                         *
*    VERIFICATION DU CALCUL DE G VIA EN 2D AVEC    *
*   ELEMENTS STANDARDS ET CHAMP_THETA UTILISATEUR  *
*                   EN MODE I PUR                  *
*                                                  *
****************************************************

* I - INIT DES DONNÉES CAS_TEST
* -----------------------------

OPTI 'DIME' 3 'ELEM' 'CUB8' ;
BTRA = FAUX ;

* DONNEES GEOMETRIQUES
A0 = 1. ;
L1 = A0 / 4. ;
DENS1 = L1 / 5. ;
DENS2 = A0 / L1 * DENS1 ;
DENS3 = 0.5*(DENS1 + DENS2) ;
DENS DENS1 ;

* PROPRIÉTÉS MATÉRIAU
MYOU = 2E11 ;
POI = 0.3 ;
KAPPA = 3-(4*POI) ;
MU = MYOU/(2*(1+POI)) ;


* II - MAILLAGE
* -------------

* CREATION DES VOLUMES
P0 = 0. 0. 0. ;
D1 = DROI (P0 MOIN (L1 0. 0.)) (P0 PLUS (L1 0. 0.)) ;
S1 = D1 TRAN (0. L1 0.) ;
CONT1 = CONT S1 ;
CONT1 = DIFF CONT1 D1 ;
CONT2 = CONT1 HOMO (A0 / L1) P0 ;
S2 = CONT1 REGL 'DFIN' DENS2 CONT2 ;
S0 = S1 ET S2 ;
V1 = S0 VOLU 'DINI' DENS3 'DFIN' DENS3 'TRAN' (0. 0. (2.*A0)) ;
V2 = V1 SYME 'PLAN' P0 (1. 0. 0.) (0. 0. 1.) ;
V0 = (V1 COUL 'BLEU') ET (V2 COUL 'ROUG') ;
* TRAC 'CACH' V0 ;FIN ;

* FUSION
PRE1 = 1.E-10 ;
PSUP = V1 POIN 'PLAN' P0 (1. 0. 0.) (0. 0. 1.) ;
PSUP = (COOR 1 PSUP) POIN 'EGSUPE' (0. - PRE1) ;
PINF = V2 POIN 'PLAN' P0 (1. 0. 0.)  (0. 0. 1.) ;
PINF = (COOR 1 PINF) POIN 'EGSUPE' (0. - PRE1) ;
ELIM PSUP PINF PRE1 ;

* FRONT
ENV0 = ENVE V0 ;
PMIL = ENV0 POIN 'PLAN' P0 (1. 0. 0.) (0. 0. 1.) ;
LEVS = ENV0 ELEM 'APPUYE' PMIL ;
* ON FAIT ATTENTION D'ORIENTER CORRECTEMENT LES LEVRES
LVSUP = LEVS ELEM 'APPUYE' (CHAN 'POI1' V1) ;
LVSUP = ORIE LVSUP (0. 1. 0.) ;
LVINF = LEVS ELEM 'APPUYE' (CHAN 'POI1' V2) ;
LVINF = ORIE LVINF (0. (-1.) 0.) ;
FRON1 = INTE (CHAN 'LIGNE' LVSUP) (CHAN 'LIGNE' LVINF) ;


* III - MODELE ET MATERIAU
* ------------------------

MOD1 = MODE V0 'MECANIQUE' 'ELASTIQUE' ;
MAT1 = MATE MOD1 'YOUN' MYOU 'NU' POI ;


* IV - CONDITIONS AUX LIMITES POUR LES 3 MODES
* --------------------------------------------

* IV.1 - EFFORTS DONNES
* *********************

* PREPARATIFS
EEXT = DIFF ENV0 LEVS ;
VEXT = V0 ELEM 'APPUYE' 'LARGEMENT' ENV0 ;
MOD2 = REDU MOD1 VEXT ;
X Y Z = COOR VEXT ;
X Y = (CHAN 'CHAM' X MOD2 'STRESSES') (CHAN 'CHAM' Y MOD2 'STRESSES') ;
TETA = CHAN (ATG Y (X + 1.E-30)) 'TYPE' 'SCALAIRE' ;
RAY1 = (((X)**2) + ((Y)**2))**0.5 ;
RAY1 = CHAN RAY1 'TYPE' 'SCALAIRE' ;
PREF = 1. / ((2*PI*RAY1)**0.5) ;
COS05 = COS (TETA/2) ;
SIN05 = SIN (TETA/2) ;
COS15 = COS (3*TETA/2) ;
SIN15 = SIN (3*TETA/2) ;
SIG0 = ZERO MOD2 'CONTRAIN' ;

* MODE I :
SXX = PREF*(COS05*(1.-(SIN05*SIN15))) ;
SXY = PREF*(COS05*SIN05*COS15) ;
SYY = PREF*(COS05*(1.+(SIN05*SIN15))) ;
SZZ = POI * (SXX + SYY) ;
SIG1 = SIG0 + (NOMC 'SMXX' SXX) + (NOMC 'SMXY' SXY) + (NOMC 'SMYY' SYY)
                + (NOMC 'SMZZ' SZZ) ;
TCHAR = REDU (BSIG MOD2 SIG1) EEXT ;


* IV.2 - DEPLACEMENT ANALYTIQUE
* *****************************

* PREPARATIFS
X Y Z = COOR V0 ;
LSUP = V1 POIN 'PLAN' P0 (1. 0. 0.) (0. 0. 1.) ;
LSUP = (COOR 1 LSUP) POIN 'EGINFE' (0. - PRE1) ;
LINF = V2 POIN 'PLAN' P0 (1. 0. 0.) (0. 0. 1.) ;
LINF = (COOR 1 LINF) POIN 'EGINFE' (0. - PRE1) ;
Y = Y + ((COOR 2 LSUP) + 1.E-30) + ((COOR 2 LINF) - 1.E-30) ;
TETA = ATG Y (X + 1.E-30) ;
RAY1 = (((X)**2) + ((Y)**2))**0.5 ;
PREF = (RAY1/(2*PI))**0.5 ;
COS05 = COS (TETA/2) ;
SIN05 = SIN (TETA/2) ;
XI0 = CHAN 'CHPO' MOD1 (ZERO MOD1 'DEPLACEM') ;

* MODE I :
XIX = PREF/(2.*MU)*(COS05*(KAPPA - 1. + (2.*(SIN05**2)))) ;
XIY = PREF/(2.*MU)*(SIN05*(KAPPA + 1. - (2.*(COS05**2)))) ;
USOL = XI0 + (NOMC 'UX' XIX) + (NOMC 'UY' XIY) ;


* V - APPEL A G_THETA
* -------------------

* G_THETA AVEC OPTION 'COUCHE'
SUPTAB = TABL ;
SUPTAB.'OBJECTIF' = MOT 'J' ;
SUPTAB.'FRONT_FISSURE' = FRON1 ; 
SUPTAB.'MODELE' = MOD1 ;
SUPTAB.'LEVRE_SUPERIEURE' = LVSUP ;
SUPTAB.'LEVRE_INFERIEURE' = LVINF ;
SUPTAB.'CARACTERISTIQUES' = MAT1 ;
SUPTAB.'SOLUTION_RESO' = USOL ;
SUPTAB.'CHARGEMENTS_MECANIQUES' = TCHAR ;
SUPTAB2 = COPI SUPTAB ;
SUPTAB.'COUCHE' = 4 ;
G_THETA SUPTAB ;
GNUM1 = SUPTAB.'CHPO_RESULTATS' ;
P1 = FRON1 POIN 1 ;
GVAL1 = EXTR (EXTR GNUM1 'VALE' P1) 1 ;

* G_THETA AVEC LE CHAMP THETA CALCULE AVEC LE NOMBRE DE COUCHES MAIS FOURNI
* DANS L'INDICE 'CHAMP_THETA'
SUPTAB2.'CHAMP_THETA' = SUPTAB.'CHAMP_THETA' ;
SUPTAB3 = COPI SUPTAB2 ;
G_THETA SUPTAB2 ;
GNUM2 = SUPTAB2.'CHPO_RESULTATS' ;
GVAL2 = EXTR (EXTR GNUM2 'VALE' P1) 1 ;

* G_THETA AVEC CHAMP THETA COMPLETEMENT CREE PAR L'UTILISATEUR 
* DANS L'INDICE 'CHAMP_THETA'
PSI PHI = PSIP V0 LVSUP 'DEUX' FRON1 ;
RHO = (((NOMC 'SCAL' PSI)**2.) + ((NOMC 'SCAL' PHI)**2.))**0.5 ;
R1 = L1 ;
R2 = 1.05*R1 ;
FUNC = BORN ((RHO - R2) / (R1 - R2)) 'COMPRIS' 0. 1. ;
THETA = MANU 'CHPO' V0 2 'UX' 1. 'UY' 0. ;
THETA = FUNC * THETA ;
SUPTAB3.'CHAMP_THETA' = THETA ;
G_THETA SUPTAB3 ;
GNUM3 = SUPTAB3.'CHPO_RESULTATS' ;
GVAL3 = EXTR (EXTR GNUM3 'VALE' P1) 1 ;

* RESULTAT ANALYTIQUE
GANA = (1. - (POI**2.)) / MYOU ;

* CALCUL DES ERREURS
SAUT 'LIGNE' ;
OPTI 'ECHO' 0 ;
MESS 'RESULTAT AVEC ''COUCHE'' :' GVAL1/50 ;
MESS 'RESULTAT AVEC ''CHAMP_THETA'' DE VALEUR IDENTIQUE :' GVAL2 ;
MESS 'RESULTAT AVEC ''CHAMP_THETA'' DIFFERENT :' GVAL3/50 ;
MESS 'RESULTAT ANALYTIQUE :' GANA/50 ;

* ERREURS
GNUM = PROG GVAL1 GVAL2 GVAL3 ;
ERR1 = GNUM - GANA ;
ERR1 = (ABS ERR1) / GANA ;
CRI1 = 1.E-2 ;
SI ((MAXI ERR1) >EG CRI1) ;
    MESS 'ERREUR : L''ERREUR SUR G DEPASSE LE CRITERE' ;
    ERRE 5 ;
FINSI ;

FIN ;