mur exposé au soleil

[ABN84]

bonsoir,
comme indiqué sur le schema ci-joint:



j'ai un mur (le grant rectangle) sur lequel est fixée une source de chaleur (le petit rectangle). ce mur est à l'equateur et dirigé est-ouest.
je voudrais determiner la temperature en chaque point de la surface de ce mur lorsque le soleil le parcoure d'est en ouest (fonction de l'angle pas du temps).
j'ai donc 5 donaines:
- contact mur/source: T1=T° cte
- contact mur/sol: T2 non cte mais invariante au cours du temps
- face est du mur: T3 variant au cours du temps
- face superieure du mur: T4 variant au cours du temps
- le reste de la face west du mur: T variant au cours du temps

merci

PS: les approximations ne me derangent pas tant que les valeurs restent acceptables et logiques
-----

[ABN84]

je remonte avant qu'il ne tombe aux oubliettes

[ABN84]

bonsoir,
ce que je cherche en realité c'est la T° en chaque point du mur (dans le volume).
j'ai ici deux programmes matlab.
le premier fonctionne, mais avec les defauts suivants:
- les conditions des 5 domaines sont choisis arbitrairement (sauf celui du radiateur qui est une donnée)
- meme nombre de pas selon hauteur et epaisseur(50) alors que la hauteur est de 2 et l'epaisseur de 0.5
- matlab inversant la representation graphique par rapport à la matrice, j'ai mon radiateur en haut et le sommet du mur en bas.

function temperature(l,h,t,max,r)
l=6;
h=2;
e=0.5;
r=50;
max=60;
lambdab=1.8;
lambdaa=0.024;
lambdat=0.75;
lambda1;
lambda2;
rho=2300;
c=800;
n=1/(r-1);
k=ones(r,r) * 273;
sigma=5.67*10^(-8);
%t0=25°C, t1=50°C
t0=298;
t1=323;
p0=t0^4*sigma;
s=l*e;
emis=0.9;
omega=ones(r,r) * 298;
d1=ones(20,1) * t1;
omega(r-19:r,r)=d1;
ps=s*emis*400*0.6/1.41;
t2=((p0+ps)/sigma)^0.25;
d2=ones(1,r) * t2;
omega(1,:)=d2;
d3=ones(30,1) * 298;
d3(1,1)=t2;
d3(30,1)=t1;
for i=2:1:29
    d3(i,1)=d3(i,1)+(-4*d3(i,1)+d3(i-1,1)+d3(i+1,1))/(2+n^2*rho*c/lambdab);
end
omega(1:30,r)=d3;
d4=ones(r,1) * t2;
omega(:,1)=d4;
d5=ones(1,r) * 298;
d5(1,1)=t2;
d5(1,50)=t1;
for j=2:1:r-1
    d5(1,j)=d5(1,j)+(-4*d5(1,j)+d5(1,j-1)+d5(1,j+1))/(2+n^2*rho*c/lambdab);
end
omega(r,:)=d5;

for i=2:1:r-1
    for j=2:1:r-1
        omega(i,j)=omega(i,j)+(-4*omega(i,j)+omega(i-1,j)+omega(i,j-1)+omega(i+1,j)+omega(i,j+1))/(4+n^2*rho*c/lambdab);
    end
end
omegac=omega-k;
contourf(omega,50);

dans le deuxieme j'ai essayé d'inverser les matrices. j'ai choisi 50 pas pour l'epaisseur et 200 pour la hauteur. les domaines d2, d3, et d4 varient en fonction de l'incliaison du soleil "a". en premier lieu je faisant comme si les surfaces etaient dissociées les unes des autres, j'obtenais donc une T° cte sur la surface verticale d2(1,200), horizontale d3(50,1), et verticale (1,180). en deuxieme lieu, je faisais intervenir la conduction de la facon suivant:
si inclinaison<45, la T° s'ecoule du bord de d2 vers d3,
si entre 45 et 90 => de d3 vers d2
si entre 90 et 135 => de d3 vers d4
si >135 => de d4 vers d3

tout ça c'est beau, mais le programme ne fonctionne pas
si qqun peut y jetter un coup d'oeuil:

max=60;
l=6;
h=2;
e=0.5;
n1=1/49;
n2=1/199;
n3=1/179;
lambdab=1.8;
lambdaa=0.024;
lambdat=0.75;
rho=2300;
c=800;
k=ones(50,200) * 273.13;
sigma=5.67*10^(-8);t0=298.13;
E1=1368;
%t1=323.13;
p0=t0^4*sigma;
s=l*e;
emis=0.9;

omega=ones(50,200) * 298.13;
d1=ones(30,1)*0;%base
d1(1,1)=298.13;
d1(30,1)=323.13;
%d2=ones(300,1)*298.15;%west
%d3=ones(1,30)*298.15;%haut
d4=ones(1,180)*298.15;%est
for m=1:1:max
        for i=2:1:49
            d1(i,1)=d1(i,1)+(-4*d1(i,1)+d1(i-1,1)+d1(i+1,1))/(2+n1*rho*c/lambdab);
        end
end
d5=ones(1,20)*50;%est/radiateur
for a=0:1:180
    if a<90
        t2=(((298.13^4)/2)+E1/(sigma*cos(a)))^0.25;
        d2=ones(1,200)*t2;%west
        t3=(((298.13^4)/2)+E1/(sigma*cos(90-a)))^0.25;
        d3=ones(50,1)*t3;%haut
        if a<45
           d3(1,1)=t2;
           for m=1:1:max
                for i=2:1:49
                   d3(i,1)=d3(i,1)+(-4*d3(i,1)+d3(i-1,1)+d3(i+1,1))/(2+n1*rho*c/lambdab);
                end
           end
        else           
           d2(1,1)=t3;
           for m=1:1:max
                for j=2:1:179
                   d2(1,j)=d2(1,j)+(-4*d2(1,j)+d2(1,j-1)+d2(1,j+1))/(2+n2*rho*c/lambdab);
                end
           end           
        end
    else
        t4=(((298.13^4)/2)+E1/(sigma*cos(180-a)))^0.25;
        d4=ones(1,180)*t4;%west
        t3=(((298.13^4)/2)+E1/(sigma*cos(a-90)))^0.25;
        d3=ones(50,1)*t3;%haut
        if a<135
           d4(0,1)=t3;
           d4(180,1)=50;
           for m=1:1:max
                for j=2:1:179
                   d4(1,j)=d4(1,j)+(-4*d4(1,j)+d4(1,j-1)+d4(1,j+1))/(2+n3*rho*c/lambdab);
                end
           end    
        else
           d3(50,1)=t4;
           for m=1:1:max
                for i=2:1:49
                   d3(i,1)=d3(i,1)+(-4*d3(i,1)+d3(i-1,1)+d3(i+1,1))/(2+n1*rho*c/lambdab);
                end
           end 
        end
    end
    d6=ones(300,1)*0;
    d6(1,180)=d4;
    d6(181:200,1)=d5;
    omega(50,:)=d6;
    omega(:,200)=d1;
    omega(1,:)=d2;
    omega(:,1)=d3;
    for j=2:1:199
        for i=2:1:49
            omega(i,j)=omega(i,j)+(-4*omega(i,j)+omega(i-1,j)+omega(i,j-1)+omega(i+1,j)+omega(i,j+1))/(4+n^2*rho*c/lambdab);
        end
    end    
    omegac=omega-k;
    %plot(omegac);
end

merci