SHARE
TWEET

Untitled

a guest Mar 10th, 2019 89 Never
Not a member of Pastebin yet? Sign Up, it unlocks many cool features!
  1. clc
  2. clear
  3. xdel(winsid())
  4.  
  5. function thetaIMenosMeio = calculaThetaRecuado(Q, pos)
  6.     if((pos == 2) | ((Q(pos)-Q(pos-1)) == 0)) then
  7.         thetaIMenosMeio = return 0;
  8.     end
  9.     thetaIMenosMeio = return (Q(pos-1)-Q(pos-2))/(Q(pos)-Q(pos-1));
  10. end
  11.  
  12. function thetaIMaisMeio = calculaThetaAvancado(Q, pos)
  13.     if((Q(pos+1)-Q(pos)) == 0) then
  14.         thetaIMaisMeio = return 0;
  15.     end
  16.     thetaIMaisMeio = return (Q(pos)-Q(pos-1))/(Q(pos+1)-Q(pos));
  17. end
  18.  
  19. function calcPsi = calculaPsiSuperbee(theta)
  20.     calcPsi = max(0, min(1, 2*theta), min(2, theta));
  21. end
  22.  
  23. function calcPsi = calculaPsiMc(theta)
  24.     calcPsi = max(0, min(((1+theta)/2), 2, 2*theta));
  25. end
  26.  
  27. function calcPsi = calculaPsiVanleer(theta)
  28.     calcPsi = (theta + abs(theta)) / (1 + abs(theta));
  29. end
  30.  
  31. function executaMetodo(calculaPsi)
  32.     I = 2000;
  33.     L = 2;
  34.     uBarra = 0.2;
  35.     deltaT = 0.005;
  36.     tempoTotal = 4;
  37.     c0 = 0.4;
  38.     C = 0.7;
  39.     deltaX = L/I;
  40.     courant = (uBarra * deltaT) / (deltaX);
  41.     disp(courant);
  42.    
  43.     qAntigo = zeros(1, I + 2);
  44.     qNovo = zeros(1, I + 2);
  45.    
  46.     i = 2;
  47.     while(i <= I+1) then
  48.         qAntigo(i) = c0;
  49.         i = i + 1;
  50.     end
  51.    
  52.     qAntigo(1) = C;
  53.     qAntigo(I+2) = qAntigo(I-1);
  54.    
  55.     tempoAtual = 0;
  56.     while(tempoAtual < tempoTotal) then
  57.         tempoAtual = tempoAtual+deltaT;
  58.         i = 2;
  59.         while(i <= I + 1) then
  60.             thetaIMenosMeio = calculaThetaRecuado(qAntigo, i);
  61.             thetaIMaisMeio = calculaThetaAvancado(qAntigo, i);
  62.             qNovo(i) = (qAntigo(i) - courant*(qAntigo(i)-qAntigo(i-1)))-0.5*courant*(1-courant)*(calculaPsi(thetaIMaisMeio)*(qAntigo(i+1)-qAntigo(i))-calculaPsi(thetaIMenosMeio)*(qAntigo(i)-                      qAntigo(i-1)));
  63.             i = i + 1;
  64.         end
  65.         qAntigo(I+2) = qAntigo(I+1);
  66.         i = 2;
  67.         while(i <= I+1) then
  68.             qAntigo(i) = qNovo(i);
  69.             i = i + 1;
  70.         end
  71.     end
  72.    
  73.     x = linspace(0, I, I);
  74.     qNovo = qNovo(1, 2:I+1);
  75.     figure();
  76.     plot2d(x, qNovo, rect=[0, 0, I, 1]);
  77.     xlabel('Número de células');
  78.     ylabel('Concentração do traçador');
  79. endfunction
  80.  
  81. executaMetodo(calculaPsiSuperbee);
  82.  
  83. executaMetodo(calculaPsiMc);
  84.  
  85. executaMetodo(calculaPsiVanleer);
RAW Paste Data
We use cookies for various purposes including analytics. By continuing to use Pastebin, you agree to our use of cookies as described in the Cookies Policy. OK, I Understand
 
Top