Risposta in frequenza su base modale

Questo articolo descrive come realizzare in Code_Aster un’analisi di risposta in frequenza su base modale. L’esempio è tratto da uno studio disponibile su internet e realizzato in Nastran. Si tratta di un elemento piastra sottoposto a due carichi: una pressione uniforme ed una forza nodale sfasata di 45 gradi. Per permettere un confronto tra i risultati verranno usati i dati dell’esempio originale (sistema imperiale). Per una descrizione più dettagliata della procedura è possibile riferirsi alla documentazione di formazione “Linear transient and harmonic analysis“.

Elemento piastra, sono rappresentati i vincoli e i carichi applicati.

Una volta creata la mesh e i gruppi necessari, è possibile procedere con la definizione del problema in Aster Study, anche in modo totalmente grafico.

Lettura e visualizzazione della mesh

La prima parte dello studio (lettura della mesh, definizione ed assegnamento dei materiali, caratteristiche degli elementi finiti, modello) procede come descritto in un post precedente. Anche la definizione dei carichi e vincoli viene fatta in modo analogo.

Prima di poter procedere con il calcolo modale occorre occuparsi della costruzione delle matrici (rigidezza, massa, smorzamento) e dei vettori relativi alle forze eccitanti (che verranno usati nell’analisi di risposta forzata). Questo viene fatto con il comando ASSEMBLAGE.

Assemblaggio delle matrici e dei vettori per l’analisi modale e la risposta forzata

Il calcolo modale si effettua con il comando CALC_MODES. Utilizziamo solo le matrici di rigidezza e di massa: in questo caso lo smorzamento verrà aggiunto direttamente durante l’analisi di risposta.

Calcolo modale tramite CALC_MODES

Prima di poter procedere con la risposta forzata, dobbiamo proiettare matrici e vettori sulla base modale, con il comando PROJ_BASE.

Proiezione sulla base modale tramite PROJ_BASE

Occorre anche definire una lista di frequenze desiderate con il comando DEFI_LIST_REEL (potremmo anche richiedere dei raffinamenti automatici con il comando DEFI_LIST_FREQ).

La lista delle frequenze su cui effettuare l’analisi di risposta

Possiamo eseguire l’analisi di risposta forzata tramite il comando DYNA_VIBRA a cui forniamo una lista delle frequenze desiderate, i due carichi dinamici (compreso lo sfasamento in gradi), le matrici, lo smorzamento modale. Indichiamo che desideriamo fare un calcolo armonico (TYPE_CALCUL=’HARM’) su base spettrale (BASE_CALCUL=’GENE’).

L’analisi di risposta in frequenza tramite DYNA_VIBRA

Essendo il calcolo effettuato su base modale, occorre riportare i risultati nello spazio fisico per poterli esaminare. Usiamo REST_GENE_PHYS a questo fine.

Conversione dei risultati tra base modale e base fisica.

A questo punto segue una manipolazione dei risultati che dipende dal singolo problema. In questo caso l’obiettivo è confrontare i risultati con quanto ci si aspetta in Nastran. Sempre dal documento iniziale troviamo per esempio la visualizzazione dello spostamento di un nodo dalla parte opposta all’eccitazione. Si nota un picco prima dei 200 Hz ed altri due intorno ai 700 dopo gli 800 Hz.

Risultati previsti, spostamento in DZ di un nodo dalla parte opposta rispetto all’eccitazione.

Possiamo estrarre i valori dello spostamento (complesso) di un nodo tramite il comando RECU_FONCTION, calcolarne il modulo tramite CALC_FONCTION e creare il grafico (o esportare in formato compatibile con Excel) tramite IMPR_FONCTION.

Risultato Code_Aster

Di seguito il file di comandi completo. La mesh e questo file possono essere scaricati da questo archivio.

# prob24
# FREQUENCY RESPONSE WITH PRESSURE AND POINT LOADS

DEBUT(PAR_LOT='NON', IGNORE_ALARM=('SUPERVIS_1'))
MAIL=LIRE_MAILLAGE(FORMAT='MED',INFO_MED=2,VERI_MAIL=_F(VERIF='OUI',),INFO=2);

MODMECA=AFFE_MODELE(MAILLAGE=MAIL,
                    AFFE=(
                          _F(GROUP_MA='CGVEGA_1',
                             PHENOMENE='MECANIQUE',
                             MODELISATION=('DKT',)),
                          ),
                    );

FCT00006=DEFI_FONCTION(
                       NOM_PARA='FREQ',
                       VALE = (
                               10, 0.1,
                               1000, 0.1,
                               ),
                       INTERPOL = ('LIN','LIN'),
                       PROL_GAUCHE='EXCLU',
                       PROL_DROITE='LINEAIRE',
                       );# Original id:310

LST00010=DEFI_LIST_REEL(
                        DEBUT = 20,
                        INTERVALLE = _F(JUSQU_A = 1000,
                                        NOMBRE = 490
                                        ),
                        );                    

# Material original id 1
M1=DEFI_MATERIAU(
                 ELAS=_F(
                         E=30000000,
                         NU=0.3,
                         RHO=0.00073038,
                         ),
                 );

# writing 0 composites
CHMAT=AFFE_MATERIAU(MAILLAGE=MAIL,
                    AFFE=(
                          _F(MATER=M1,GROUP_MA=('CGVEGA_1',),),
                          ),
                    );

CAEL=AFFE_CARA_ELEM(MODELE=MODMECA,
                    COQUE=(
                           _F(GROUP_MA='CGVEGA_1',
                              EPAIS=0.1,
                              VECTEUR=(0.9,0.1,0.2)),
                           ),
                    );

# ConstraintSet original id : 1
BL2=AFFE_CHAR_MECA(MODELE=MODMECA,
                   DDL_IMPO=(
                             _F(NOEUD=('N1', 'N4', 'N24', 'N35', 'N46', ),DX=0, DY=0, DZ=0, DRX=0, DRY=0, ),
                             ),
                   );

# LoadSet original id : 100
CHMEC3=AFFE_CHAR_MECA(MODELE=MODMECA,
           PRES_REP=(
                         _F(PRES= 1,GROUP_MA='CGVEGA_1',),
                      ),
                      );

# LoadSet original id : 600
CHMEC6=AFFE_CHAR_MECA(MODELE=MODMECA,
                      FORCE_NODALE=(
                                    _F(NOEUD='N21',FZ=1,),
                                    ),
                      );

# Analysis original id : 1
ASSEMBLAGE(MODELE=MODMECA,
           CHAM_MATER=CHMAT,
           CARA_ELEM=CAEL,
           CHARGE=(
                   BL2,
                   ),
           NUME_DDL=CO('NUMDDL1'),
           MATR_ASSE=(_F(OPTION='RIGI_MECA', MATRICE=CO('RIGI1'),),
                      _F(OPTION='MASS_MECA', MATRICE=CO('MASS1'),),
                      _F(OPTION='AMOR_MECA', MATRICE=CO('AMOR1'),),
                      ),
          VECT_ASSE=(
                      _F(OPTION='CHAR_MECA',VECTEUR=CO('FX1_5'),CHARGE=(CHMEC3,)),
                      _F(OPTION='CHAR_MECA',VECTEUR=CO('FX1_6'),CHARGE=(CHMEC6,)),
                      ),
           );

MODES1=CALC_MODES(MATR_RIGI=RIGI1,
                       MATR_MASS=MASS1,
                       SOLVEUR_MODAL=_F(METHODE='TRI_DIAG'),
                       OPTION='BANDE',
                       CALC_FREQ=_F(
                                    FREQ=(10, 2000),
                                    ),
                       VERI_MODE=_F(STOP_ERREUR='NON',),
                       SOLVEUR=_F(METHODE='MUMPS',
                                  RENUM='PORD',
                                  NPREC=8,
                                  ),
                       );

PROJ_BASE(BASE=MODES1,
          MATR_ASSE_GENE=(_F(MATRICE=CO('MASSG1'),
                             MATR_ASSE=MASS1,),
                          _F(MATRICE=CO('RIGIG1'),
                             MATR_ASSE=RIGI1,),
                          _F(MATRICE=CO('AMORG1'),
                             MATR_ASSE=AMOR1,),
                          ),
          VECT_ASSE_GENE=(
                          _F(VECTEUR=CO('VG1_5'),
                             VECT_ASSE=FX1_5,),
                          _F(VECTEUR=CO('VG1_6'),
                             VECT_ASSE=FX1_6,),
                         ),
          );

GENE1 = DYNA_VIBRA(
                   TYPE_CALCUL='HARM',
                   BASE_CALCUL='GENE',
                   MATR_MASS  = MASSG1,
                   MATR_RIGI  = RIGIG1,
                   MATR_AMOR  = AMORG1,
                   AMOR_MODAL = _F(AMOR_REDUIT = 0.03,),
                   LIST_FREQ  = LST00010,
                   EXCIT      = (
                                 _F(
                                    VECT_ASSE_GENE = VG1_5,
                                    FONC_MULT = FCT00006,
                                    PHAS_DEG = 0,),
                                 _F(
                                    VECT_ASSE_GENE = VG1_6,
                                    FONC_MULT = FCT00006,
                                    PHAS_DEG = -45,),
                                 ),
                   );

RESU1 = REST_GENE_PHYS(RESU_GENE = GENE1,
                       TOUT_ORDRE = 'OUI',
                       NOM_CHAM = 'DEPL',
                       );

FODZ55=RECU_FONCTION(RESULTAT = RESU1,
                      NOM_CMP = 'DZ',
                      NOEUD = 'N55',
                      NOM_CHAM = 'DEPL');                    
                      

FODZ55M=CALC_FONCTION(EXTRACTION=_F(FONCTION=FODZ55,PARTIE='MODULE'),);

IMPR_FONCTION(FORMAT='XMGRACE',
              UNITE=27,
              ECHELLE_Y='LOG',
              PILOTE = 'PNG',
              COURBE=(_F(FONCTION=FODZ55M,COULEUR=3,MARQUEUR=0,)))                          

IMPR_RESU(FORMAT='MED',
          UNITE=80,
          RESU=(_F(RESULTAT=RESU1,PARTIE='REEL',NOM_CHAM="DEPL",NOM_CHAM_MED="DEPLM"),
                _F(RESULTAT=RESU1,PARTIE='IMAG',NOM_CHAM="DEPL",NOM_CHAM_MED="DEPLP")));

FIN(FORMAT_HDF='OUI',);

Lascia un commento

Il tuo indirizzo email non sarà pubblicato. I campi obbligatori sono contrassegnati *

Questo sito usa Akismet per ridurre lo spam. Scopri come i tuoi dati vengono elaborati.