Calcolo del primo modo proprio di vibrare

Si vuole ora calcolare il primo modo proprio di vibrare della membrana e la corrispondente lunghezza d'onda. Il problema da risolvere risulta il seguente:

$$\frac{\partial^2 u}{\partial x^2} + \frac{\partial^2 u}{\partial y^2} + \kappa^2 u = 0$$ (22)

con le condizioni al contorno definite dalla (2). Nella (22) $\kappa$ è funzione dello spostamento $u$ ed è inversamente proporzionale alla lunghezza d'onda. Dalla teoria del calcolo variazionale si determina che le lunghezze d'onda ammissibili per la membrana in esame corrispondono ai punti di stazionarietà, cioè i punti con gradiente nullo, del funzionale:

$$\kappa^2 \left( u \right) = \frac{\displaystyle{\int \!\!\! \int_... ... \!\!\! \int_{\cal S}} \;\displaystyle{ \frac{1}{2}} u^2 \left( x,y \right) dS}$$ (23)

Fra tutti i punti di stazionarietà posseduti dal funzionale (23), noi siamo interessati al primo, corrispondente al valore minimo di $\kappa$ associato al modo di vibrare della membrana con lunghezza d'onda massima.

Il funzionale (23) si calcola numericamente usando la soluzione approssimata $\hat{u}$ della (6). Si osservi che $\hat{u}$ puø essere riscritta come il seguente prodotto scalare:

$$\hat{u} = \mbox{\boldmath$\xi$} ^T \bu$$ (24)

dove $\xi$ è il vettore delle funzioni forma $\left[ \xi_1, \xi_2, \ldots, \xi_n \right]^T$. Utilizzando la (24) il numeratore della (23) risulta:

$$\int \!\!\! \int_{\cal S} \; \frac{1}{2} \left[ \left( \frac{\par... ...l y} \frac{\partial \mbox{\boldmath$\xi$}^T}{\partial y} \right] dS \right) \bu$$ (25)

Le funzioni forma $\xi_k$ hanno supporto locale e pertanto la matrice definita dall'espressione integrale che compare nella (25) va calcolata mediante la consueta procedura di assemblaggio dei contributi derivanti dai termini locali. L'elemento in riga $i$ e colonna $j$ di tale matrice risulta pertanto:

$$\sum_e \int \!\!\! \int_{{\cal S}^{(e)}} \; \left[ \frac{\partial... ...\xi_i} {\partial y} \frac{\partial \xi_j}{\partial y} \right] dS^{(e)} = h_{ij}$$ (26)

cioè proprio l'elemento in riga $i$ e colonna $j$ della matrice di rigidezza $H$ (vedi equazione (10)). Il numeratore del funzionale (23) approssimato mediante la (6) si può dunque scrivere come:

$$\int \!\!\! \int_{\cal S} \; \frac{1}{2} \left[ \left( \frac{\par... ...artial \hat{u}}{\partial y} \right)^2 \right] dS = \frac{1}{2} \; \bu^T \bH \bu$$ (27)

Procedendo in modo analogo per il denominatore si ottiene:

$$\int \!\!\! \int_{\cal S} \: \frac{1}{2} \hat{u}^2 dS = \frac{1}{... ...\!\! \int_{\cal S} \mbox{\boldmath$\xi$} \mbox{\boldmath$\xi$}^T dS \right) \bu$$ (28)

L'espressione integrale definisce una matrice $P$ il cui generico elemento $p_{ij}$, ottenuto mediante l'assemblaggio dei contributi locali, risulta:

$$p_{ij} = \sum_e p_{ij}^{(e)} = \sum_e \int \!\!\! \int_{{\cal S}^{(e)}} \xi_i \xi_j dS^{(e)}$$ (29)

Pertanto, il denominatore si scrive come:

$$\int \!\!\! \int_{\cal S} \: \frac{1}{2} \hat{u}^2 dS = \frac{1}{2} \bu^T \bP \bu$$ (30)

e la matrice $P$ è denominata matrice di capacità della membrana. Ricordando che:

$$\int \!\!\! \int_{{\cal S}^{(e)}} \xi_i \xi_j dS^{(e)} = \left\{ ... ... \displaystyle{\frac{\Delta}{12}} \;\;\;\;\;\; i \neq j \end{array} \right.$$ (31)

si ottiene la matrice di capacità locale su un generico elemento triangolare:

$$P^{(e)} = \frac{\Delta}{12} \left[ \begin{array}{ccc} 2 & 1 & 1 1 & 2 & 1 1 & 1 & 2 \end{array} \right]$$ (32)

Dalla (29) si può osservare che $P$ è una matrice simmetrica ed è facile verificare che è anche definita positiva.

La procedura per assemblare le $P^{(e)}$ nella $P$ è identica a quella implementata per assemblare le $H^{(e)}$ nella $H$. Si noti anche che la topologia di $P$ sarà identica a quella di $H$, vale a dire che i vettori JA ed IA non cambiano.

L'espressione numerica del funzionale (23), ottenuta discretizzando numeratore e denominatore secondo la (27) e la (30), è dunque quella del quoziente di Rayleigh nella forma generalizzata:

$$\lambda \left( \hat{u} \right) = \kappa^2 \left( \hat{u} \right) = \frac{\bu^T \bH \bu}{\bu^T \bP \bu}$$ (33)

I vettori che rendono stazionaria la (33) rappresentano i modi di vibrare della membrana e $\lambda$ l'inverso del quadrato della corrispondente lunghezza d'onda. La quantità cercata $\lambda_n$ è pertanto l'autovalore minimo del problema generalizzato:

$$H \bu = \lambda P \bu$$ (34)

ed il modo proprio di vibrare il corrispondente autovettore $\bu_n$. Tali quantità possono essere efficientemente determinate minimizzando il quoziente di Rayleigh generalizzato (33) col metodo del GCM. Anche in questo caso, affinché il problema sia determinato, è necessario imporre le condizioni al contorno (2) (membrana fissa al bordo) sia su $P$ che su $H$. Si osservi che, tuttavia, l'imposizione delle condizioni al contorno mediante la procedura indicata in precedenza modifica lo spettro del problema generalizzato (34) in funzione dei valori $R_{max}$ prescelti. Per evitare di operare sugli ultimi autovalori, e quindi di modificare il risultato del problema, è consigliabile utilizzare nella matrice $P$ un valore di $R_{max}$ di qualche ordine di grandezza inferiore a quello adottato in $H$.

La convergenza viene raggiunta quando il residuo relativo risulta inferiore ad una prefissata tolleranza. Si può, ad esempio, adottare la medesima tolleranza utilizzata per la soluzione del sistema lineare o leggermente superiore nel caso di una convergenza più difficoltosa.