Il Bar dell'Ingegneria

Alessia LT

« Older   Newer »
 
  Share  
.
  1.     +1   -1
     
    .
    Avatar

    Advanced Member

    Group
    Administrator
    Posts
    8,169
    Reputation
    +295

    Status
    Offline
    Apro questo nuovo topic, che riguarda il tema delle funi-membrane-tensostrutture, per poter costituire contenitore del programma AlessiaLT e di quanto gentilmente l'ing E. Casagrande ci vorrà illustrare in merito.

    Naturalmente questo topic puo' anche essere sede di discussione.


    rQXT6fo


    L'ing. E.Casagrande ha gia condiviso l'eseguibile del suo programma in altro topic "Dynamic Relaxation BarbaTrick" e comunque potete scaricarlo dal link che segue:

    Download AlessiaLT
     
    Top
    .
  2. ecasagrande
        +1   -1
     
    .

    User deleted


    Ringrazio Afazio per l'apertura di questo post.

    Iniziamo con il dire che ALESSIA LT è un puro solutore con un visualizzatore di risultati in formato testo. Questo per sottolineare che, per chi mi ha contattato, non si vede il classico form con puntatoti etc., ma esso è un motore privo di carozzeria!
    La sigla LT sta ad indicare, invece, una versione Light. Questa versione è stata posta di tipo LT perchè il DRM ha molti parametri di stabilità, i quali ovviamente influenzano la stabilità di calcolo. Tra questi si citano:
    - l'intervallo di tempo;
    - il coefficiente di smorzamento.
    L'intervallo di tempo è stato assunto costante e pari a 1.01; tale valore risulta stimato in base a quanto specificato nel libro a pag. 75 (il paragrafo è abbastanza impegnativo); si potrebbero introdurre altri valori ma prima di far ciò si deve essere consci del calcolo in sé.
    La tolleranza è stata ssunta pari a 0.0001 con iterazione massima pari a 10000. In effetti questi valori portrebbero anch'essi essere variati, ma aumentare il numero di iterazioni potrebbe influenzare il processo numerico.
    CODICE
    INTEGER, PARAMETER             :: nmax=10000
         REAL*8,    PARAMETER           :: toll=0.0001, dt=1.01, massa=1000000000000.0

    L'algoritmo segue il flusso logico presentato a pag. 81 che successivamente presenterò.
     
    Top
    .
  3.     +1   -1
     
    .
    Avatar

    Advanced Member

    Group
    Administrator
    Posts
    8,169
    Reputation
    +295

    Status
    Offline
    Immagino che quando riporti il numero di pagina ti riferisci al tuo libro.
    Ho comprato il tuo libro, ma l'ho a casa. Potrei scansionare le pagine a cui ti riferisci ed eventualmente postarle o postarne uno stralcio?
     
    Top
    .
  4. ecasagrande
        +1   -1
     
    .

    User deleted


    Allora, per questioni burocratiche non posso riportare nulla (non è cattiveria ma spero che potete capirmi!). Facciamo in questo modo. Se siete d'accordo, cerco di spiegare a grandi linee in modo da non appesantire la trattazione che è già di per sé molto pesante.....
     
    Top
    .
  5.     +1   -1
     
    .
    Avatar

    Advanced Member

    Group
    Member
    Posts
    2,942
    Reputation
    +187

    Status
    Offline
    Un ri-benvenuto a Casagrande. Che in questo secondo incontro con noi vedo più 'fecondo' di interesse.

    Avendoci 'sbattuto' la testa non possiamo che essere contenti dei tuoi interventi.
    Visto che ci dici che Alessia LT è solamente 'motore' e non 'carrozzeria', e visto che vedo il file di input come semplice file di testo, volevo chiederti: anche il file di output è un file di testo?

    Mi piego meglio.
    Il 'nostro' visualizzatore si vorrebbe farlo diventare prima o poi 'processore' di qualcosa.
    Al momento ci interessava validare con afazio quanto avevamo pensato circa la sola visualizzazione di un modello 3D sia in vista assonometrica che prospettica. Il risultato ci ha sorpreso positivamente noi per primi.
    Proprio perchè inizialmente si volevano testare le sole opzioni di visualizzazione, per far prima, ho semplicemente creato un file di testo da cui il 'visualizzatore' legge i dati che gli servono (coordinate dei nodi, funi, ecc.). Nella fattispecie i nodi, le funi, ecc. sono quelli di una 'semplice' struttura che afazio ha risolto in excel con una procedura semi-automatica.
    Nel seguito io per primo voglio aggiungere al visualizzatore una 'corposa' serie di opzioni di editing. Ma il tempo è sempre tiranno....

    Detto questo, poichè attualmente il visualizzatore non fa altro che leggere file di testo, mi basta sapere il 'formato' del file di output sfornato da Alessia LT per poter far visualizzare gli esiti del calcolo con Alessia.
     
    Top
    .
  6.     +1   -1
     
    .
    Avatar

    Advanced Member

    Group
    Administrator
    Posts
    8,169
    Reputation
    +295

    Status
    Offline
    CITAZIONE (ecasagrande @ 10/2/2014, 11:49) 
    Allora, per questioni burocratiche non posso riportare nulla (non è cattiveria ma spero che potete capirmi!). Facciamo in questo modo. Se siete d'accordo, cerco di spiegare a grandi linee in modo da non appesantire la trattazione che è già di per sé molto pesante.....

    Capisco benissimo (diritti dell'editore) non fartene un cruccio. Risolviamo come hai delineato.
     
    Top
    .
  7. ecasagrande
        +1   -1
     
    .

    User deleted


    CITAZIONE
    Un ri-benvenuto a Casagrande. Che in questo secondo incontro con noi vedo più 'fecondo' di interesse

    Lo ero anche prima forse non lo davo da vedere in modo marcato!

    CITAZIONE
    Detto questo, poichè attualmente il visualizzatore non fa altro che leggere file di testo, mi basta sapere il 'formato' del file di output sfornato da Alessia LT per poter far visualizzare gli esiti del calcolo con Alessia.

    Si Alessia LT genera un file output tipo quello allegato nel file zip (lo allego nuovamente) in cui ci sono le informazioni con i risultati.

    Allora passiamo al trattazione. Il passo temporale dt dovrebbe essere, secondo studi specifici, minore di 2/(periodo della struttura) e pertanto funzione di K e M (rigidezza e massa). In base a valutazioni che ho effettuato in diverse strutture il dt ottimale si aggira intorno a 1.0 e 1.1. Attenzione che se si inserisce anche un piccolo discostamento il metodo potrebbe diventare instabile! Ad esempio per una struttura semplice l'assumere 1.1 al posto che 1.0 porta ad ottenere una soluzione in 183 iterazioni anzichè in 213. Pare una nullità ma provate a pensare quando i nodi sono migliaia.
    Un altro metodo potrebbe essere quello di Qiang il quale aggiorna il dt ad ogni passo in base al periodo di vibrazione della struttura. Sistema ottimale che prevede due calcoli distinti ad ogni iterazione: calcolo dinamico con valutazione del primo modo e calcolo non lineare. Ovviamente i tempi sono lunghissimi in quanto la rigidezza aumenta sempre più.

    Per lo smorzamento, invece, la cosa è assai più complicata. E' per questo che è stato introdotto il kinetic damping che ha il compito di annullare il coefficiente di smorzamento. C'è un altro parametro però molto importante ed è la massa. Infatti il DRM è fondato sul comportamento dinamico della struttura. La massa può essere assunta proporzionale ad un coefficiente, oppure alla rigidezza. Di solito scelgo il metodo di Gerschgorin che in pratica afferma che la massa deve essere maggiore o uguale a 1/4 dt^2 Somma(Kij).
    File Allegato
    output.dat
    (Number of downloads: 94)

     
    Top
    .
  8. ecasagrande
        +1   +1   -1
     
    .

    User deleted


    Scelti i parametri iniziali possiamo passare all'implementazione del metodo. Si precisa che la scelta dei parametri potrebbe essere molteplice grazie all'inserimento di subroutine specifiche richiamate all'occorenza.
    A questo punto allego un diagramma dello schema logico per l'implementazione del DRM con la variante "kinetic damping". Acquisiti i dati iniziali si passa al calcolo non lineare. Una cosa va detta per i dati di input. Sinceramente per la "pulizia" dell'algoritmo sarebbe preferibile l'inserimento dei nodi vincolati per primi; questo permette una più agevole lettura (a parere mio) dell'intero processo di calcolo.
    Come si nota dal diagramma, il calcolo non lineare viene inizializzato attraverso un ciclo while:
    CODICE
    !inizio del calcolo non lineare
         iter=1
         err=1000.0
         DO WHILE (err.GE.toll)
            !inizializza la rigidezza concentrata
            DO i=1,nj
               rigid(i,1)=0.0
            ENDDO

    A questo punto si deve calcolare la matrice tangente e calcolare il residuo che vedremo prossimamente.
    Attached Image
    Diagramma

     
    Top
    .
  9.     +1   -1
     
    .
    Avatar

    Advanced Member

    Group
    Administrator
    Posts
    8,169
    Reputation
    +295

    Status
    Offline
    Fin qui tutto chiaro anche perchè avevo sviscerato la questione in autonomia.
    L'unica cosa che non capisco è la fase che hai denominato "calcolo matrice di rigidezza".

    In cosa consiste nel dettaglio e come viene determinata?
     
    Top
    .
  10. ecasagrande
        +1   -1
     
    .

    User deleted


    Con il blocco "matrice di rgidezza" ho indicato il processo per calcolare il contributo di ogni singolo elemento che compone il sistema. Il DRM si basa su l'approccio del residuo determinato come differenza tra forze esterne e forze interne (già ampiamente visto). Da tale assunzione ne deriva che il termine di confronto deve essere fatto su ogni singolo nodo; e qui viene il bello.
    La rigidezza quindi, non viene calcolata come i classici metodi FEM, ovvero con la creazione della matrice globale, ma viene concentrata "lumped" al nodo. In questo modo la rigidezza in un generico nodo risulta pari al contributo delle singole aste afferenti al nodo in questione.
    Il calcolo viene eseguito preliminarmente determinando il contributo tangente KG il quale deve essere sommato successivamente al contributo elastico KE. Effettuando quindi una ricorsione per ogni nodo e per ogni asta si arriva alla definizione della rigidezza al nodo con cui proseguire il calcolo.

    Effettuato questo processo si parte con il calcolo delle forze interne.....
     
    Top
    .
  11.     +1   -1
     
    .
    Avatar

    Advanced Member

    Group
    Member
    Posts
    2,942
    Reputation
    +187

    Status
    Offline
    Casagrande, in qualche prova fatta da lui, alla fine Afazio era giunto alla conclusione che il DRM non potrebbe essere utilizzato nella fase iniziale di ricerca di forma della struttura, ma solamente per le fasi successive in cui, con forma già trovata, si vanno direttamente a caricare i nodi.

    Per la fase iniziale si dovrebbero utilizzare altre metodologie di calcolo, come quello che fa uso del Metodo della densità di forza e della matrice di connettività.

    Tu confermi questa sua 'conclusione'?
     
    Top
    .
  12.     +1   -1
     
    .
    Avatar

    Advanced Member

    Group
    Administrator
    Posts
    8,169
    Reputation
    +295

    Status
    Offline
    CITAZIONE (ecasagrande @ 12/2/2014, 15:30) 
    Con il blocco "matrice di rgidezza" ho indicato il processo per calcolare il contributo di ogni singolo elemento che compone il sistema. Il DRM si basa su l'approccio del residuo determinato come differenza tra forze esterne e forze interne (già ampiamente visto). Da tale assunzione ne deriva che il termine di confronto deve essere fatto su ogni singolo nodo; e qui viene il bello.
    La rigidezza quindi, non viene calcolata come i classici metodi FEM, ovvero con la creazione della matrice globale, ma viene concentrata "lumped" al nodo. In questo modo la rigidezza in un generico nodo risulta pari al contributo delle singole aste afferenti al nodo in questione.
    Il calcolo viene eseguito preliminarmente determinando il contributo tangente KG il quale deve essere sommato successivamente al contributo elastico KE. Effettuando quindi una ricorsione per ogni nodo e per ogni asta si arriva alla definizione della rigidezza al nodo con cui proseguire il calcolo.

    Effettuato questo processo si parte con il calcolo delle forze interne.....

    Ok, Chiaro.
    Ma nel caso di ricerca della forma iniziale come ci si comporta dato che in genere partiamo da una configurazione arbitraria che puo' dar luogo a rigidezze nulle secondo una direzione?
    Mi riferisco, per esempio, al caso in cui partiamo da una maglia di nodi posti tutti ad una prefissata quota, tranne ovviamente i nodi vincolati che hanno una quota di progetto fissata.

    In questo caso considerando un nodo intermedio e non adiacente ad un nodo vincolato, tutte le funi ad esso convergenti saranno orizzontali e non ammettendo nessuna componente verticale determinano rigidezza nulla nella direzione z.
     
    Top
    .
  13. ecasagrande
        +1   -1
     
    .

    User deleted


    Cerco di rispondere sia a Zax2013 che ad Afazio.
    Allora, concordo con quanto discusso da entrambi sul problema del form-finding utilizzando il DRM (tant'è vero che nella maggior parte dei casi utilizzo altri algoritmi per la ricerca dello stato "zero"). Basti pensare all'implementazione del DRM; infatti esso si basa sul concetto "dinamico" ovvero ho una struttura in un certo stato, muovo le masse e definisco una condizione di quiete associata ad uno smorzamento. Risulta ovvio che se ho uno stato anomalo, in cui la struttura è in equilibrio orizzontale risulta alquanto difficile la ricerca di forma.
    Ma qui voglio sottolineare una cosa che ho evidenziato anche nel libro. Molti ipotizzano un peso proprio della fune nullo; questo in realtà non è assolutamente vero in quanto come ben sapete dalla semplice risoluzione dell'equazione della catenaria: peso nullo consegue un cavo rettilineo con pretensione nulla! Se alla struttura invece si carica il peso concentrato della fune le cose ovviamente cambiano.
    C'è un'altra cosa da dire. Per sistemi uniformi esistono anche accortezze particolari. Mi spiego meglio. Per effettuare la ricerca di forma di una rete uniforme quadrangolare è bene osservare quanto segue:
    - per mantenere l'equidistanza si deve considerare una rigidezza molto elevata; questa sembra una banalità ma è fondamentale in quanto il DRM è un metodo dinamico facilmente instabile;
    - molti utilizzatori effettuano il form-finding agendo sul sistema di vincolamento. In pratica si sostituisce al vincolo un elemento con bassa rigidezza e lo si collega ad un altro vincolo. Questi elementi possono essere caratterizzati da una tensione costante oppure essere rimossi controllando la forma.
    Sono procedure complesse che non sempre sono attuabili. Prova ne sia che in forme complesse le aste vengono raggrupate per controllare meglio il calcolo.
     
    Top
    .
  14. ecasagrande
        +1   -1
     
    .

    User deleted


    Nell'ambito del DRM è bene porre l'attenzione su un'altra questione. Nella variante con salto dell'energia cinetica, il DRM deve calcolare il valore dell'energia e qualora esso sia di "picco", la velocità viene inizializzata con valore nullo. Nel calcolo però non è ben stabilito se il valore ottenuto è effettivamente di picco, ma esso evidenzia solamente è stato sorpassato il valore più alto. Ebbene considerando che i risultati sono funzione di tale parametro, succede molte volte che la deformata in genere, si discosti un pò da altri metodi (rigidezze). Per ovviare a questo, molti utilizzando un passo temporale molto piccolo in modo tale da centrale il massimo; questo sistema, ovviamente, non è ottimale sia perchè rende il calcolo instabile (ricordate dt 1.00 - 1.01) e perchè si aumentano in modo sconsiderato i tempi di calcolo.
    Per ovviare a tale inconveniente, sono disponibili degli algoritmi specifici proprio per riceercare il picco. Tali metodi si basano sulla ricerca di un punto massimo su una funzione quadratica tra due intervalli di tempo operando una sorta di "retromarcia" controllata.
    Si riporta di seguito il codice che implementa tale correzione.
    CODICE
    DO i=1,gdl
               keold=keold+0.5*(mt(1,i)*v1(i,1)**2)
               kenew=kenew+0.5*(mt(1,i)*n(i,1)**2)
            ENDDO
            IF (kenew.LT.keold) THEN
               DO i=1,gdl
                  cor(i,1)=-(3*dt*n(i,1)/2)+(dt**2*r(i,1))/(2*m(i,1))
                  disp(i,1)=disp(i,1)+cor(i,1)
               ENDDO
     
    Top
    .
  15.     +1   -1
     
    .
    Avatar

    Advanced Member

    Group
    Administrator
    Posts
    8,169
    Reputation
    +295

    Status
    Offline
    CITAZIONE (ecasagrande @ 15/2/2014, 09:46) 
    Nell'ambito del DRM è bene porre l'attenzione su un'altra questione. Nella variante con salto dell'energia cinetica, il DRM deve calcolare il valore dell'energia e qualora esso sia di "picco", la velocità viene inizializzata con valore nullo. Nel calcolo però non è ben stabilito se il valore ottenuto è effettivamente di picco, ma esso evidenzia solamente è stato sorpassato il valore più alto. Ebbene considerando che i risultati sono funzione di tale parametro, succede molte volte che la deformata in genere, si discosti un pò da altri metodi (rigidezze). Per ovviare a questo, molti utilizzando un passo temporale molto piccolo in modo tale da centrale il massimo; questo sistema, ovviamente, non è ottimale sia perchè rende il calcolo instabile (ricordate dt 1.00 - 1.01) e perchè si aumentano in modo sconsiderato i tempi di calcolo.
    Per ovviare a tale inconveniente, sono disponibili degli algoritmi specifici proprio per riceercare il picco. Tali metodi si basano sulla ricerca di un punto massimo su una funzione quadratica tra due intervalli di tempo operando una sorta di "retromarcia" controllata.
    Si riporta di seguito il codice che implementa tale correzione.
    CODICE
    DO i=1,gdl
               keold=keold+0.5*(mt(1,i)*v1(i,1)**2)
               kenew=kenew+0.5*(mt(1,i)*n(i,1)**2)
            ENDDO
            IF (kenew.LT.keold) THEN
               DO i=1,gdl
                  cor(i,1)=-(3*dt*n(i,1)/2)+(dt**2*r(i,1))/(2*m(i,1))
                  disp(i,1)=disp(i,1)+cor(i,1)
               ENDDO

    ehehe "retromarcia controllata"

    Nel mio codice in VBA avevo scritto proprio per controllare la retromarcia, il seguente codice:

    CODICE
    ' riporto le coordinate dei punti a quelle dello step precedente
           For i = 1 To nTot
               X(i) = X(i) - 3 * Vx(i) * DT / 2 + 0.5 * DT ^ 2 * Rx(i) / m
               Y(i) = Y(i) - 3 * Vy(i) * DT / 2 + 0.5 * DT ^ 2 * Ry(i) / m
               Z(i) = Z(i) - 3 * Vz(i) * DT / 2 + 0.5 * DT ^ 2 * Rz(i) / m
           Next


    In cui pero le componenti di velocità venivano determinate attraverso il seguente:

    CODICE
    Sub CalcolaVelocita(start As Boolean)
    Dim i As Integer

    For i = 1 To nTot
       If start Then
           Vx(i) = (1 - vinc(i)) * 0.5 * DT * Rx(i) / m
           Vy(i) = (1 - vinc(i)) * 0.5 * DT * Ry(i) / m
           Vz(i) = (1 - vinc(i)) * 0.5 * DT * Rz(i) / m
       Else
           Vx(i) = (1 - vinc(i)) * (Vx(i) + DT * Rx(i) / m)
           Vy(i) = (1 - vinc(i)) * (Vy(i) + DT * Ry(i) / m)
           Vz(i) = (1 - vinc(i)) * (Vz(i) + DT * Rz(i) / m)
       End If
    Next
    End Sub

    Il parametro "start" passato è True quando siamo nelle condizioni di azzeramento dell'energia cinetica.
     
    Top
    .
18 replies since 10/2/2014, 10:39   1617 views
  Share  
.