Accueil | Cast3M

* fichier : iss3D_xyz.dgibi
************************************************************************
************************************************************************
*
************************************************************************
*           REPONSE SISMIQUE DU SOL EN ABSENCE DE STRUCTURE  
************************************************************************
*
* DESCRIPTION DU PROBLEME
* ------------------------
* IL S'AGIT D'UN PROBLEME D'INTERACTION SOL-STRUCTURE. 
* EN ABSENCE DE STRUCTURE, IL N'Y PAS D'INTERACTION. ON DOIT RETROUVER
* L'ACCELEROGRAMME IMPOSE A LA SURFACE DU SOL VIA LE PROCESSUS DE 
* DECONVOLUTION (PROCEDURES DECONV OU DECONV3D) ET CONVOLUTION 
* (PROCEDURE DYNAMIC OU PASAPAS).
* LE CYLINDRE DANS CE JEU DE DONNEES REPRESENTE LE SOL PROCHE TANDIS 
* QUE LE SOL LOINTAIN QUI S'ETEND VERS L'INFINI EST REPRESENTE PAR UNE 
* FRONTIERE ABSORBANTE COMPOSEE D'AMORTISSEURS VISQUEUX.
************************************************************************
*-----------------------------------------------------------------------
*              PROCEDURE DE CALCUL DE LA DERIVEE PREMIERE
*-----------------------------------------------------------------------
'DEBPROC' DERIV1 EV1 ;
     X     = 'EXTR' EV1 'ABSC'  ;
     Y     = 'EXTR' EV1 'ORDO'  ;
     NPOIN = 'DIME' X ;
     H     = ( 'EXTR' X 2 ) - ( 'EXTR' X 1 ) ;
 
     X1    = ( 'PROG' 0.0 0.0 ) 'ET' Y ;
     X3    = Y 'ET' ( 'PROG' 0.0 0.0 ) ;
     DX    = ( X3 - X1 ) / ( H * 2.0 ) ;
     Y1    = 'ENLE' ( 'ENLE' DX ( NPOIN + 2 ) ) 1 ;
     L1    = ( 'PROG' 0.0 ) 'ET' ( 'PROG' ( NPOIN - 2 ) * 1.0 )
                            'ET' ( 'PROG' 0.0 ) ;
     Y2    = Y1 * L1 ;
     EV2   = 'EVOL' 'MANU'  X Y2  ;
'FINPROC' EV2 ;
*-----------------------------------------------------------------------
 
*-----------------------------------------------------------------------
OPTI TRAC PSC ;
OPTI DIME 3 ELEM CU20 ;
 
R1  = 2.0 / (PI**0.5) ;
H0  = 0.0 ;
NT  = 2  ;
ANG = 90. / NT ;
NE  = 2 ; 
 
P0 = 0   0  H0 ;
P1 = R1  0  H0 ;
P2 = TOUR P1 (0 0 0) (0 0 1) ANG ;
*-----------------------------------------------------------------------
*                      SOL SANS STRUCTURE (CU20, PR15)
*-----------------------------------------------------------------------
L01B = DROI NE P0 P1 ;
L12B = DROI 1  P1 P2 ;
L20B = DROI NE P2 P0 ;
 
L1B = L01B ET L12B ET L20B ;
ELIM L1B 0.001 ;
 
S1B  = COUL VERT (SURF PLAN L1B) ;
SQ1B = S1B ;
 
I = 1 ;
REPE B1B (NT - 1) ;
S2B  = TOUR S1B (0 0 0) (0 0 1) ANG ;
SQ1B = SQ1B ET S2B ;
S1B  = S2B ;
I = I + 1 ;
FIN B1B ;
 
NR = 3 ;
R2 = NR*R1;
 
P4 = R2  0  H0 ;
P5 = TOUR P4 (0 0 0) (0 0 1) ANG ;
 
L14 = DROI ((NR - 1)*NE) P1 P4 ;
L45 = DROI         1     P4 P5 ;
L52 = DROI ((NR - 1)*NE) P5 P2 ;
L21 = DROI         1     P2 P1 ;
 
SS1  = DALL L14 L45 L52 L21 PLAN ;
SSQ1 = SS1 ;
 
I = 1 ;
REPE B2 (NT - 1) ;
SS2  = TOUR SS1 (0 0 0) (0 0 1) ANG ;
SSQ1 = SSQ1 ET SS2 ;
SS1 = SS2 ;
I = I + 1 ;
FIN B2 ;
 
SSQ1 = SSQ1 ET SQ1B ;
ELIM SSQ1 0.001 ;
TRAC SSQ1 QUAL NCLK ;
 
SSQ2 = SYME SSQ1 DROIT P0 (0 1 H0) ;
SSQ3 = SYME SSQ1 POINT P0 ;
SSQ4 = SYME SSQ1 DROIT P0 P1 ;
 
SSOL = (SSQ1 ET SSQ2 ET SSQ3 ET SSQ4) COUL VERT ;
 
SOL1 = VOLU TRANS ((NR/2)*NE) SSOL (0 0 (NR*R1/-2)) ;
ELIM SOL1 0.001 ;
TRAC SOL1 CACH NCLK ;
 
SOL2 = SOL1 PLUS (0 0 (NR*R1/-2)) ;
 
TITR 'MAILLAGE DU SOL' ;
SOL = SOL1 ET SOL2 ;
ELIM SOL 0.001 ;
 
TITR 'MAILLAGE DU SOL, NB ELEMENT =' (NBEL SOL) ;
TRAC SOL CACH NCLK ;
 
ENV1 = ENVE SOL1 ;
ENV2 = ENVE SOL2 ;
 
HS1 = NR*R1/-2 ;
HS2 = NR*R1/-1 ;
 
S11 = ENV1'POIN' 'PLAN' (0 0 H0) (1 0 H0) (0 1 H0) 0.001     ;
S11 = ENV1'ELEM' 'APPUYE' 'STRICTEMENT' S11                  ;
 
S12 = ENV1'POIN' 'PLAN' (0 0 HS1) (1 0 HS1) (0 1 HS1) 0.001  ;
S12 = ENV1'ELEM' 'APPUYE' 'STRICTEMENT' S12                  ;
 
S21 = ENV2 'POIN' 'PLAN' (0 0 HS1) (1 0 HS1) (0 1 HS1) 0.001 ;
S21 = ENV2 'ELEM' 'APPUYE' 'STRICTEMENT' S21                 ;
 
S22 = ENV2 'POIN' 'PLAN' (0 0 HS2) (1 0 HS2) (0 1 HS2) 0.001 ;
S22 = ENV2 'ELEM' 'APPUYE' 'STRICTEMENT' S22                 ;
 
BOR1 = DIFF (DIFF ENV1 S11) S12 ;
BOR2 = DIFF (DIFF ENV2 S21) S22 ;
 
FOND_S = S22 ;
 
TRAC (BOR1 ET BOR2 ET FOND_S) NCLK ;
*-----------------------------------------------------------------------
*                              MODELE SOL
*-----------------------------------------------------------------------
NU1  = 1./3. ;
G1   = 100.E7 ;
E1   = 2*(1 + NU1)*G1 ;
MESS 'E =' E1  ; 
 
MOD_S  = MODE SOL MECANIQUE ELASTIQUE ISOTROPE CU20 PR15 ;
MAT_S  = MATE MOD_S YOUN E1 NU NU1 RHO 2000 ;
MAS_S  = MASS MOD_S MAT_S ;
RIG_S  = RIGI MOD_S MAT_S ;
*-----------------------------------------------------------------------
*                         AMORTISSEMENT RAYLEIGH
*-----------------------------------------------------------------------
F1 = 5.0 ;
F2 = 25.0 ;
KSI_S   = 0.05 ;
ALPHA   = 4.0 * PI * F1 * F2 / ( F1 + F2 ) ;
BETA    = 1.0 / ( PI * ( F1 + F2 ) ) ;
AMO_S   = KSI_S  * ((ALPHA * MAS_S ) ET (BETA * RIG_S )) ;
*-----------------------------------------------------------------------
*              SIGNAL SISMIQUE (SELON SPECTRE PS92 SITE S1)
*-----------------------------------------------------------------------
LFR  = PROG 0.1 0.5 1.0 2.5 5.  33. 50. ; 
LSP  = PROG 0.1 0.5 1.0 2.5 2.5 1.0 1.0 ; 
SP0  = EVOL MANU 'FREQ(HZ)' LFR 'ACCE(M/S*S)' LSP ;
 
TAB = TABLE ;
TAB.'MOTIT' = 'SPECTRE PS92 S1 ' ;
TAB.'SEISME'= TABLE ;
TAB.'SEISME'.'SPECTRE'  =  SP0 ;
TAB.'SEISME'.'AMORT'    =  0.05 ;
TAB.'SEISME'.'TYPSP'    = 'ACCE' ;
TAB.'SIGNAL' = TABLE ;
TAB.'SIGNAL'.'ENVE'     = 'PLATLIN' ;
TAB.'SIGNAL'.'NP'       =  8 ;
TAB.'SIGNAL'.'DUREE'    =  2.56 ;
TAB.'SIGNAL'.'TDEBUT'   =  0.75 ;
TAB.'SIGNAL'.'TFIN'     =  1.5 ;
TAB.'NBITER'            =  5 ;
TAB.'NBSIGN'            =  3 ;
TAB.'NALEAT'            =  3 ;
TAB.'FRCOUP'            =  49.9 ;
TAB.'OPTSORT'           = 'SPECTRE' ;
 
TABSIG = SIGNSYNT FABR TAB ;
 
*   MOYENNE ZERO, PLAGE INITIALE ZERO
 
LT0 = EXTR TABSIG.1 ABSC ;
LA0 = EXTR TABSIG.1 ORDO ;
DT  = (EXTR LT0 2) - (EXTR LT0 1) ;
NP0 = DIME LT0 ;
LT1 = PROG 0. PAS DT 2.3 ;
LA1 = IPOL LT1 LT0 LA0 ;
NP1 = DIME LT1 ; 
LA1 = LA1 - (PROG NP1 * ((SOMM LA1) / NP1)) ;
LA1 = (PROG 20 * 0.) ET LA1 ;
LT1 = PROG 0. PAS DT NPAS (19 + NP1) ;
ACCX = EVOL MANU 'TEMPS(S)' LT1 'ACCE(M/S*S)' LA1 ;
 
LT0 = EXTR TABSIG.2 ABSC ;
LA0 = EXTR TABSIG.2 ORDO ;
DT  = (EXTR LT0 2) - (EXTR LT0 1) ;
NP0 = DIME LT0 ;
LT1 = PROG 0. PAS DT 2.3 ;
LA1 = IPOL LT1 LT0 LA0 ;
NP1 = DIME LT1 ; 
LA1 = LA1 - (PROG NP1 * ((SOMM LA1) / NP1)) ;
LA1 = (PROG 20 * 0.) ET LA1 ;
LT1 = PROG 0. PAS DT NPAS (19 + NP1) ;
ACCY = EVOL MANU 'TEMPS(S)' LT1 'ACCE(M/S*S)' LA1 ;
 
LT0 = EXTR TABSIG.3 ABSC ;
LA0 = EXTR TABSIG.3 ORDO ;
DT  = (EXTR LT0 2) - (EXTR LT0 1) ;
NP0 = DIME LT0 ;
LT1 = PROG 0. PAS DT 2.3 ;
LA1 = IPOL LT1 LT0 LA0 ;
NP1 = DIME LT1 ; 
LA1 = LA1 - (PROG NP1 * ((SOMM LA1) / NP1)) ;
LA1 = (PROG 20 * 0.) ET LA1 ;
LT1 = PROG 0. PAS DT NPAS (19 + NP1) ;
ACCZ = (EVOL MANU 'TEMPS(S)' LT1 'ACCE(M/S*S)' LA1) *2./3. ;
 
DESS ACCX MIMA NCLK ;
DESS ACCY MIMA NCLK ;
DESS ACCZ MIMA NCLK ;
*-----------------------------------------------------------------------
*                 DECONVOLUTION ET FRONTIERE ABSORBANTE
*-----------------------------------------------------------------------
TAB = TABLE ;
 
TAB.1 = TABLE ;
TAB.1 .'FRONTIERE'        = BOR1  ;
TAB.1 .'MASSE_VOLUMIQUE'  = 2000. ;
TAB.1 .'POISSON'          = 1./3. ;
TAB.1 .'YOUNG'            = E1    ;
TAB.1 .'AMORTISSEMENT'    = 0.05  ;
 
TAB.2 = TABLE ;
TAB.2 .'FRONTIERE'        = BOR2  ;
TAB.2 .'MASSE_VOLUMIQUE'  = 2000. ;
TAB.2 .'POISSON'          = 1./3. ;
TAB.2 .'YOUNG'            = E1    ;
TAB.2 .'AMORTISSEMENT'    = 0.05  ;
 
P_GAM  = TABLE ;
P_GAM . 1 = P4 ;
 
TABX = DECONV3D TAB P0 FOND_S MOD_S MAT_S 'UX' 
                ACCX 5. 25. 50.  P_GAM ;
TABY = DECONV3D TAB P0 FOND_S MOD_S MAT_S 'UY' 
                ACCY 5. 25. 50.  P_GAM ;
TABZ = DECONV3D TAB P0 FOND_S MOD_S MAT_S 'UZ' 
                ACCZ 5. 25. 50.  P_GAM ;
*-----------------------------------------------------------------------
*                  TABLE POUR LA PROCEDURE 'DYNAMIC'
*-----------------------------------------------------------------------
CH_DEPI = MANU CHPO 3 SOL UX 0.0 UY 0.0 UZ 0.0 ;
CH_VITI = MANU CHPO 3 SOL UX 0.0 UY 0.0 UZ 0.0 ;
 
TAB_DYN = TABLE ;
TAB_DYN.'DEPL' = CH_DEPI ;
TAB_DYN.'VITE' = CH_VITI ;
TAB_DYN.'RIGI' = RIG_S ;
TAB_DYN.'MASS' = MAS_S ;
TAB_DYN.'AMOR' = AMO_S ET TABX.'AMOR' ;
TAB_DYN.'CHAR' = TABX.'CHAR' ET TABY.'CHAR' ET TABZ.'CHAR' ;
TAB_DYN.'FREQ' = TABX.'FCDYN' ;
DT             = 1 * TABX.'PAS' ;
NB_PAS         = ENTI (2.51 / DT) ;
TAB_DYN.'DEBU' = 0.0 ;
TAB_DYN.'INST' = PROG 0.0 PAS DT NPAS (NB_PAS - 1) ;
 
TDYNA = DYNAMIC TAB_DYN ;
*-----------------------------------------------------------------------
*                          POST-TRAITEMENT
*-----------------------------------------------------------------------
POI = TABLE ;
POI.1 = P0  ;
POI.2 = P4  ;
 
NP = DIME POI ;
 
L_TEM = PROG NB_PAS * 0.0 ;
L_VX = TABLE ;
L_VY = TABLE ;
L_VZ = TABLE ;
 
I = 1 ;
REPE B1 NP ;
L_VX.I  = PROG NB_PAS * 0.0 ;
L_VY.I  = PROG NB_PAS * 0.0 ;
L_VZ.I  = PROG NB_PAS * 0.0 ;
I = I + 1 ;
FIN B1 ;
 
I_TEM = 0 ;
REPETER B2 NB_PAS ;
  I_TEM = I_TEM + 1           ;
  TEM_I = TDYNA.I_TEM.'TEMP'   ;
  REMP L_TEM I_TEM TEM_I      ;
  CHV_I = TDYNA.I_TEM.'VITE'   ;
 
  I = 1 ;
    REPE B2_1 NP ;
    REMP L_VX.I I_TEM (EXTR CHV_I POI.I UX) ;
    REMP L_VY.I I_TEM (EXTR CHV_I POI.I UY) ;
    REMP L_VZ.I I_TEM (EXTR CHV_I POI.I UZ) ;
  I = I + 1 ;
  FIN B2_1 ;
 
FIN B2 ;
 
*   REPONSE EN ACCELERATION 
 
LT1 = EXTR ACCX ABSC ;
EV_VITX = TABLE ;
EV_VITY = TABLE ;
EV_VITZ = TABLE ;
 
EV_ACCX = TABLE ;
EV_ACCY = TABLE ;
EV_ACCZ = TABLE ;
 
I = 1 ;
REPE B3 NP ;
EV_VITX.I = EVOL MANU L_TEM L_VX.I ;
EV_VITY.I = EVOL MANU L_TEM L_VY.I ;
EV_VITZ.I = EVOL MANU L_TEM L_VZ.I ;
 
EV_ACCX.I = DERIV1 EV_VITX.I ;
EV_ACCY.I = DERIV1 EV_VITY.I ;
EV_ACCZ.I = DERIV1 EV_VITZ.I ;
 
LT2 = EXTR EV_ACCX.I ABSC ;
LA2 = EXTR EV_ACCX.I ORDO ;
LA2 = IPOL LT1 LT2 LA2 ;
EV_ACCX.I = EVOL MANU 'TEMPS(S)' LT1 'ACCE(M/S*S)' LA2 ;
 
LT2 = EXTR EV_ACCY.I ABSC ;
LA2 = EXTR EV_ACCY.I ORDO ;
LA2 = IPOL LT1 LT2 LA2 ;
EV_ACCY.I = EVOL MANU 'TEMPS(S)' LT1 'ACCE(M/S*S)' LA2 ;
 
LT2 = EXTR EV_ACCZ.I ABSC ;
LA2 = EXTR EV_ACCZ.I ORDO ;
LA2 = IPOL LT1 LT2 LA2 ;
EV_ACCZ.I = EVOL MANU 'TEMPS(S)' LT1 'ACCE(M/S*S)' LA2 ;
 
I = I + 1 ;
FIN B3 ;
 
*   COMPARAISON AVEC L'ACCELEROGRAMME IMPOSE
 
SI (NEG GRAPH 'N');
DESS ((COUL ROUG EV_ACCX.1) ET (COUL VERT EV_ACCX.2) ET ACCX) NCLK 
TITR 'ACCELEROGRAMMES X IMPOSE (BLEU) ET CALCULES (ROUGE ET VERT)' MIMA ;
FINSI ;
 
DESS ((COUL ROUG EV_ACCY.1) ET (COUL VERT EV_ACCY.2) ET ACCY) NCLK
TITR 'ACCELEROGRAMMES Y IMPOSE (BLEU) ET CALCULES (ROUGE ET VERT)' MIMA ;
FINSI ;
 
DESS ((COUL ROUG EV_ACCZ.1) ET (COUL VERT EV_ACCZ.2) ET ACCZ) NCLK
TITR 'ACCELEROGRAMMES Z IMPOSE (BLEU) ET CALCULES (ROUGE ET VERT)' MIMA ;
FINSI ;
 
ERR1 = MAXI ABS (EXTR (EV_ACCX.1 - ACCX) ORDO) ;
ERR2 = MAXI ABS (EXTR (EV_ACCX.2 - ACCX) ORDO) ;
ERR3 = MAXI ABS (EXTR (EV_ACCY.1 - ACCY) ORDO) ;
ERR4 = MAXI ABS (EXTR (EV_ACCY.2 - ACCY) ORDO) ;
ERR5 = MAXI ABS (EXTR (EV_ACCZ.1 - ACCZ) ORDO) ;
ERR6 = MAXI ABS (EXTR (EV_ACCZ.2 - ACCZ) ORDO) ;
 
SI ((MAXI (PROG ERR1 ERR2 ERR3 ERR4 ERR5 ERR6)) > 0.01);
ERRE 5;
FINSI;
 
FIN ;