Discussion :

[emagis12]

bonjour tout le monde;
Merci de me faire comprendre le programme suivant "a propos de la simulation" comment ca fonctionne; merci beaucoup :

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
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
 
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% cavity2d.m: 2D cavity flow, simulated by a LB method            
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Lattice Boltzmann sample, Matlab script
% Copyright (C) 2006-2008 Jonas Latt
% Address: Rue General Dufour 24,  1211 Geneva 4, Switzerland 
% E-mail: Jonas.Latt@cui.unige.ch
%
% Implementation of 2d cavity geometry and Zou/He boundary
% condition by Adriano Sciacovelli
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% This program is free software; you can redistribute it and/or
% modify it under the terms of the GNU General Public License
% as published by the Free Software Foundation; either version 2
% of the License, or (at your option) any later version.
% This program is distributed in the hope that it will be useful,
% but WITHOUT ANY WARRANTY; without even the implied warranty of
% MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
% GNU General Public License for more details.
% You should have received a copy of the GNU General Public 
% License along with this program; if not, write to the Free 
% Software Foundation, Inc., 51 Franklin Street, Fifth Floor,
% Boston, MA  02110-1301, USA.
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 
clear
 
% GENERAL FLOW CONSTANTS 
lx = 128; 
ly = 128; 
 
uLid  = 0.1; % horizontal lid velocity 
vLid  = 0;    % vertical lid velocity 
Re    = 100;  % Reynolds number 
nu    = uLid *lx / Re;     % kinematic viscosity 
omega = 1. / (3*nu+1./2.); % relaxation parameter 
maxT  = 40000; % total number of iterations 
tPlot = 10;    % cycles for graphical output
 
% D2Q9 LATTICE CONSTANTS 
t   = [4/9, 1/9,1/9,1/9,1/9, 1/36,1/36,1/36,1/36]; 
cx  = [ 0, 1, 0, -1, 0, 1, -1, -1, 1]; 
cy  = [ 0, 0, 1, 0, -1, 1, 1, -1, -1]; 
opp = [ 1, 4, 5, 2,  3, 8, 9,  6,  7]; 
lid = [2: (lx-1)]; 
 
[y,x] = meshgrid(1:ly,1:lx); 
obst = ones(lx,ly); 
obst(lid,2:ly) = 0; 
bbRegion = find(obst); 
 
% INITIAL CONDITION: (rho=0, u=0) ==> fIn(i) = t(i) 
fIn = reshape( t' * ones(1,lx*ly), 9, lx, ly); 
 
% MAIN LOOP (TIME CYCLES) 
for cycle = 1:maxT 
 
% MACROSCOPIC VARIABLES 
rho = sum(fIn); 
ux = reshape ( (cx * reshape(fIn,9,lx*ly)), 1,lx,ly ) ./rho; 
uy = reshape ( (cy * reshape(fIn,9,lx*ly)), 1,lx,ly ) ./rho; 
 
% MACROSCOPIC (DIRICHLET) BOUNDARY CONDITIONS 
 
ux(:,lid,ly) = uLid; %lid x - velocity 
uy(:,lid,ly) = vLid; %lid y - velocity 
rho(:,lid,ly) = 1 ./ (1+uy(:,lid,ly)) .* ( ... 
                sum(fIn([1,2,4],lid,ly)) + 2*sum(fIn([3,6,7],lid,ly)) ); 
 
% MICROSCOPIC BOUNDARY CONDITIONS: LID (Zou/He BC)
fIn(5,lid,ly) = fIn(3,lid,ly) - 2/3*rho(:,lid,ly).*uy(:,lid,ly); 
fIn(9,lid,ly) = fIn(7,lid,ly) + 1/2*(fIn(4,lid,ly)-fIn(2,lid,ly))+ ... 
                1/2*rho(:,lid,ly).*ux(:,lid,ly) - 1/6*rho(:,lid,ly).*uy(:,lid,ly); 
fIn(8,lid,ly) = fIn(6,lid,ly) + 1/2*(fIn(2,lid,ly)-fIn(4,lid,ly))- ... 
                1/2*rho(:,lid,ly).*ux(:,lid,ly) - 1/6*rho(:,lid,ly).*uy(:,lid,ly); 
 
% COLLISION STEP 
for i=1:9 
    cu = 3*(cx(i)*ux+cy(i)*uy); 
    fEq(i,:,:) = rho .* t(i) .* ... 
        ( 1 + cu + 1/2*(cu.*cu) - 3/2*(ux.^2+uy.^2) ); 
    fOut(i,:,:) = fIn(i,:,:) - omega .* (fIn(i,:,:)-fEq(i,:,:)); 
end 
 
% MICROSCOPIC BOUNDARY CONDITIONS: NO-SLIP WALLS (bounce-back)
for i=1:9 
    fOut(i,bbRegion) = fIn(opp(i),bbRegion); 
end 
 
% STREAMING STEP 
for i=1:9 
    fIn(i,:,: ) = circshift(fOut(i,:,: ), [0,cx(i),cy(i)]); 
end
 
% VISUALIZATION
if (mod(cycle,tPlot)==0)
    u = reshape(sqrt(ux.^2+uy.^2),lx,ly);
    uxx=reshape(ux,lx,ly);
    uyy=reshape(uy,lx,ly);
    u(bbRegion) = nan;
    Unorm=u(:,ly:-1:1)'./uLid;
    maxU=max(Unorm);minU=min(Unorm);
    imagesc(Unorm);
    title(['contour velocity     Nbr cycle=' num2str(cycle) ')']);
    colorbar
    axis equal off; drawnow
end
 
end

[duf42]

Bonjour,

Pourrais-tu préciser ta question? Quelles sont les lignes du code que tu ne comprends pas exactement?

Duf

[emagis12]

bjr,ca va bien?
ben je fais des etudes en simulation fluide:par exemple simuler un ecoulement de fluide ou de chaleur et dorenavant je dois ecrire des programmes matlab destines a des simulations pareilles! j'aimerais k vous m'aidez a decortiquer tout ce programme ,les differentes etapes de ce programmes et a koi ca servent chacunes de ces boucles. merci d'avance

[duf42]

On veut bien t'aider mais il va falloir que tu commences à défricher un peu tout ca de ton côté.

On ne peut pas t'expliquer l'intégralité du programme d'un coup, d'une part car celà représente un travail titanesque et d'autre part car nous n'avons pas forcément les compétences techniques pour interpréter correctement l'application du programme à la mécanique des fluides.

[ysagon]

Bonjour, ce code provient de notre site web palabos
Vous en avez juste enlevé le copyright..