Discussion :

[Antoine10]

Bonsoir,
Je suis actuellement en 3ème année de physique et je dois rendre un projet sur matlab (celui ci est codé scilab) qui est le suivant :

On dispose de 16 dipôles magnétiques qui sont répartis sur un réseau cubique (4*4) avec pour paramètre de maille a.
Autrement dit, vous prenez une grille de 4 sur 4 et en chacun des 16 points d'intersection, vous mettez un dipôle magnétique avec pour moment magnétique m.
Tous ces dipôles ont leur moment magnétique qui est orienté vers le "haut" (soit un angle initial de 0 radian) et vous mettez uniquement lemoment magnétique du 16ème dipôle orienté d'un angle pi/2 par rapport aux 15 autres. On note téta (fonction du temps) cet angle.

L'objectif est de trouver la valeur des 16 tétas en fonction du temps afin de pouvoir créer une simulation qui montrerai l'évolution du mouvement des dipôles.

J'ai déjà pas mal avancé sur le projet mais les résultats trouvés ne me plaisent pas beaucoup...
Je trouve un graphe qui montre l'évolution des 16 angles tétas en fonction du temps mais il me paraît faux.
Par exemple, lorsqu'on teste un angle de 0 radian comme condition initiale pour le 16ème dipôle (à la place de pi/2), logiquement on ne devrait pas observer de mouvement puisque les 16 dipôles sont en équilibres (enfin je pense...)

Je sais que cela peut demander un peu de temps mais si quelqu'un pouvait m'aider à trouver ce qui cloche ce serait top !
Je met mon script en dessous.
A mon avis, le problème se trouve dans la résolution de l'équation différentielle (peut être pourrait on utiliser ode ...)

La formule pour calculer le champ magnétique créé par les dipôle est une formule connue (vous pouvez me la demander si vous n'êtes pas sur de son expression). De plus, dans l'équa diff, on apraxie la dérivée seconde par une différence finie (méthode alternative qui remplace ode)

Merci beaucoup pour votre aide.
Cordialement, Antoine.

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
3
4
function [R] = fonction_distance( i, X, Y) // Calcul de la distance entre la source et l'observateur
 
R=sqrt(((X-X(i))).^2+((Y-Y(i)).^2));
R(i)=[];

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
3
4
function [R5] = fonction_invdist( i, X, Y) // Calcul de l'inverse de la distance puissance 5 entre la source et l'observateur
 
R=fonction_distance( i, X, Y);
R5=1 ./ (R.^5);

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
function [B] = champ_magnetique(teta, X, Y, a) // Calcul du champ magnétique total aux 16 points du réseau
 
i=1;
while i<=16
    R=fonction_distance( i, X, Y)
    R5= fonction_invdist( i, X, Y );
 
    T=teta;
    T(i)=[];
    C=cos(T);
    S=sin(T);
    Ax=(X'-X(i))./a;
    By=(Y'-Y(i))./a;
    Ax(i)=[];
    By(i)=[];
 
    B(i,1)= (3.*a^5.*R5*((S*Ax+C*By)*Ax-(((R.^2)./3)*diag(S))')./(a^2));
    B(i,2)= (3.*a^5.*R5*((S*Ax+C*By)*By-(((R.^2)./3)*diag(C))')./(a^2));
    i=i+1;
end

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
function [D] = fonction_diff( teta, tetam1, X, Y, T0, Tmax, dt, a)
 
NP=round(Tmax/dt); // NP nombre de points total
 
D=(zeros(NP,16)); // Taille donnée aux variables de la boucle
C=zeros(16,16);
S=zeros(16,16);
 
B=champ_magnetique( teta, X, Y, a); // calcul du 1er champ magnétique
 
for k=[3 : NP]
    D( k,: )=teta;
    C=cos(teta);
    S=sin(teta);
    tetap1=((dt/T0)^2).*((B(1:16,2)')*diag(S)-(B(1:16,1)')*diag(C))-tetam1+2.*teta;
    B=champ_magnetique(tetap1, X, Y, a);
    tetam1=teta;
    teta=tetap1;
end

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
clear all
 
u0=4*%pi*10^-7; // perméabilité magnétique du vide
 
a=1; // paramètre du réseau cubique
Y = a.*[0 0 0 0 1 1 1 1 2 2 2 2 3 3 3 3];
X = a.*[0 1 2 3 0 1 2 3 0 1 2 3 0 1 2 3]; // référentiel du réseau
N = (X+4.*Y)+1; // Numérotation de la maille
 
m=1*10^(-10); // norme du moment magnétique
J=1; // moment d'inertie
 
B0=u0*m/(4*%pi*(a^3));
T0=sqrt(J/(m*B0)); // Normalisation de l'échelle
 
Nperiode=500; // point par période
Ntot=1; // nombre de périodes
Tmax=Ntot*T0;
dt=T0/Nperiode; 
 
teta=zeros(1,15); 
teta=[teta %pi/2];
tetam1=teta; // conditions initiales
 
D=fonction_diff( teta, tetam1, X, Y, T0, Tmax, dt, a);
 
t=[0:dt:Tmax-dt]; // vecteur temps
 
plot(t,D) // tête des tetas en fontion du temps

[Michel]

Bonjour Antoine10,
je ne suis pas un spécialiste de scilab et c'est en cherchant des exemples que je suis tombé sur ton problème où tu écris
Par exemple, lorsqu'on teste un angle de 0 radian comme condition initiale pour le 16ème dipôle (à la place de pi/2), logiquement on ne devrait pas observer de mouvement puisque les 16 dipôles sont en équilibres (enfin je pense...)
Pour moi , non !
Les dipôles extérieurs n'ont pas de voisins d'un côté. Cette particularité explique la disparition progressive du magnétisme dans un barreau de fer aimanté.