Résolution numérique : EDP d'ordre 4

[erff]

Bonjour !

Voici mon problème : je veux résoudre numériquement cette EDP


Avec

EDIT : je n'arrive pas à éliminer le "n" apres le arctan : ne pas en tenir compte

L'idée avec g est de faire une fonction |signe(u0''(x))| de manière régulière, V est un paramètre qui permet d'adoucir la transition de 1 à 0 aux pts d'annulation de la dérivée seconde.


- Tout d'abord j'ai cru comprendre qu'il manquait des conditions aux bords pour que le pb admette une seule solution ? Est ce le cas ?


- Ensuite, j'ai essayé de résoudre ceci par différences finies (Euler à droite, à gauche, Crank-Nicholson) mais aucune des solution ne converge à moins de mettre un pas de temps de 0.001 et un pas de 1 sur les x...vraiment pas pratique.

Pour Euler à droite j'ai un truc du genre :
où U^n représente le vecteur colonne u(x) à t=t_n
A est une matrice ayant la diagonale principale et les 4 au dessus de remplies.

- Un autre pb est qu'avec C-N, je ne sais pas quoi remplir dans les 4 premières et 4 dernières lignes des matrices (A et B ci dessous), car même avec des conditions aux bords on n'arrrive pas à rendre inversible cette matrice

(relation du genre :



Pour résumer : voici mes questions

- Est ce normal que la solution soit presque tout le temps divergente avec ces méthodes numériques ?

- Quelle méthode numérique utiliser pour une telle EDP (moi je ne connais que les différences finies qui marchaient pas trop mal pr l'équation de la chaleur) ?


Et au passage, vous auriez un bouquin à me conseiller qui traîte de la résolution des EDP (avec un peu de théorie quand même) pour un niveau en maths de math spé (bac+2) ?


Merci beaucoup
-----

[edpiste]

Pour que le problème soit bien posé, il faut effectivement mettre des conditions de bord. Sur quel intervalle (en x) veux-tu résoudre ton EDP ? R tout entier ?
Et si tu as essayé d'implémenter une méthode numérique, tu as dû implicitement écrire des conditions de bord (quel intervalle as-tu pris en x ? quelles valeurs pour u aux extrémités de cet intervalle ? quelle donnée initiale ?)

Si j'étais toi, je commencerai par essayer de faire converger ma méthode numérique préférée dans le cas simple g(x)=1 (toutes celles que tu as citées doivent marcher dans ce cas).
Pour ma part, je ferais sans doute un schéma d'Euler implicite en temps avec résolution d'un problème elliptique à chaque étape par différences finies standard.

[erff]

Bonjour, et merci pour ta réponse.

Depuis, j'ai reprogrammé une résolution numérique par Crank-Nicholson et je me suis rendu compte que j'avais fait beaucoup d'erreur dans l'écriture de la matrice (je ne suis pas très habitué à programmer ).


- x varie sur un intervalle de longueur fixée (fini)
- Mon objectif est d'atteindre u(t=+oo,x) (ou du moins une approximation)

En gros, voilà ce que j'ai fait :

- dérivée temporelle évaluée par : [u(n+1)-u(n)]/dt (pour x fixé)

- Calcul du terme d/dx ( g(x)*d^3u/dx^3 ) qui donne une relation de récurrence entre u(j-2),...,u(j+2)) : je calcule toujours les dérivée de sorte à rester "centré", je ne sais pas si je suis clair...

- Calcul de ce terme avec Euler implicite donc la relation de récurrence est de la forme : où n représente l'instant t_n.

- Idem avec Euler explicite :

- Donc on évalue ce terme en faisant la moyenne :



- Je ne connais grand chose aux EDP donc j'ai refait ce qu'on avait vu une fois en cours pour calculer une dérivée seconde : on prenait celle calculée au temps n, celle calculée au temps n+1 et on faisait la moyenne des 2...j'ai fait exactement pareil pour une dérivation d'ordre 4.

Au final, j'obtient le système matriciel :



Où A et B sont des matrices qui ont les diagonales "de -2 à +2 (notation matlab)" remplies et qui se ressemblent beaucoup (juste la diagonale qui diffère d'un terme +2*dx^4/dt et -2*dx^4/dt).

Mes conditions aux limites sont les suivantes :
- Dérivée nulle aux bords (en x=0 et x=K)
- Dérivée troisième nulle aux bords (idem)

Pour ce faire, je résouds un systeme de taille K-4*K-4 (K=taille de U^n, pour tout n) et je calcule à part les termes u1,u2, u(K-1) et u(K)

Résultats :

- Le truc ne diverge jamais (enfin !!!), mais a un comportement bizarre si l'on prend dx^4 petit devant dt (comme si on changeait le signe à chaque itération) ce qui se comprend car dans ce cas les matrices A et B sont "presque égales".

Maintenant je me demande si ça converge vers la bonne solution, et quelles sont les conditions optimales ? car mon but est d'obtenir la solution asymptotique, peu importe les étapes intermédiaires.


Voilà voilà, j'espère avoir tout dit !

PS : Si vous avez des schémas numériques à conseiller je suis ouvert


Merci bcp.

[erff]

Je me permet d'ajouter quelques résultats de simulation :

- Est ce que cela est plausible ? J'ai bien l'impression que oui, car mon but était de créer un signal "qui reflète" les variations lentes du signal initial en conservant les points d'inflexion, ce qui semble être le cas ici.

- Est-ce normal que l'algo ne s'arrete jamais de vraiment converger c à d que si j'itère bcp bcp de fois, la solution change qud meme entre les itérations (de manière infime entre les itérations (quasiment invisible), mais en cumulant tout ça, l'effet est non négligeable) ? Peut on justifier cela avec les erreurs inhérentes au schéma numérique ?

Voilà !

Ci joints mes résultats de simulation

[A voir en vidéo sur Futura]

[edpiste]

Pas eu le temps de tout lire dans les détails. J'essaierai de te faire une réponse plus complète une prochaine fois. Je suis un peu étonné que la solution ne tende pas vers 0, lorsque t tend vers l'infini. En dehors des points où la dérivée seconde de u0 s'annule, l'équation est dissipative et devrait tout écraser.
Que donne ton schéma dans le cas où g(x)=1 ? Si la solution ne tend pas vers zéro, il y a un pb. Aussi, quel est l'intérêt de faire la moyenne entre Euler implicite et explicite ? Ca améliore la convergence ?

[erff]

Bonjour,

Mon souci est résolu, apparemment tout était lié à une instabilité numérique (pentes raides pour la fonction g). Cependant en prenant g=1 le truc convergeait vers 0 effectivement, mais la simple présence de 0 dans g empêchait le truc de diffuser. J'ai donc changé d'EDP (bah oui qud il y a pas de solution, on change le problème ! ) pour me ramener à ceci :

Finalement ceci converge vraiment bien et préserve bien les points d'inflexions.

Aussi, quel est l'intérêt de faire la moyenne entre Euler implicite et explicite ? Ca améliore la convergence ?

C'est la méthode de Cranck-Nicholson, elle est inconditionnellement stable (enfin je l'avais vue pour les EDP d'ordre 2, j'ai adapté à l'ordre 4 et ça a l'air de bien fonctionner).
J'ai aussi pris pour g une fonction qui préserve les maxima et le truc marche super bien. Cela dit le côté théorique reste encore un peu obscur pour moi.

Merci pour tes interventions.

[erff]

J'en profite pour demander quelques éclaircissements :

- Y a t-il un moyen de déterminer un ordre de grandeur du pas de temps dt à ne pas dépasser en fonction de dx pour que le schéma converge raisonnablement vers la solution ? Car cette résolution coûte cher en temps de calcul, et comme je l'ai précisé, je ne veux garder que la solution asymptotique (pour t->+oo) : les étapes intermédiaires ne m'intéressent pas.

Merci