» »

Matlab problem

Matlab problem

kitajc ::

V Matlabu sem izračunal enačbo za 2 podani vrednosti (input =...). Sedaj pa bi rad, da mi ti dve vrednosti pobira avtomatsko in izračuna končne podatke.

Se pravi:
Pr, kjer so vrednosti od 0.01; 0.05; 0.1; 0.2; 0.4; 0.6; 0,8; 1; 1.2; 1.5; 2; 3; 5; 7; 10
Tr, kjer so vrednosti od 0.30; 0.35; 0.40 ...... 4.00 (korak med njimi pa ni vseskozi enak)

Izračunati pa bi mi moral vse možne kombinacije. Zaporedje pri tem ni pomembno, ker bom rezultate vstavil v tabelo. Tako, da bom lahko na koncu pogledal - pri Tr = 0.4 in Pr = 1, je rezultat enačbe "tak".

Bi kdo to znal nastaviti?

kitajc ::

Npr. imam enačbo:

X = a + b
vrednosti za a so 1; 3; 6; 7; 10; 12;
vrednosti za b so 2; 3; 8; 17; 25; 44;

Kako lahko nastavim, da ne bom vnašal vrednosti vsakič posebej? Zelo lepo bi bilo dobiti rezultat v obliki tabele ali matrike.

chort ::

Jaz bi se tega lotil takole: kreiraš dva vektorja, prvi vsebuje vse a-je, drugi vse b-je.
Nato bi z FOR zanko šel skozi prvi vektor a[i] za vse i-je, vgnezdeno v to pa še FOR zanko ki gre skozi drugi vektor b[j] za vse j-je. Znotraj te druge zanke bi izvajal tvoj ukaz (X = a + b).

Če uporabiš indeksa s katerimi se drajsaš skozi oba zgornja vektorja (recimo jima i in j)tudi kot definicijo lokacije v matriki za izhodne podatke X[i,j], boš dobil matriko z rezultati.
Srečno.

bluefish ::

Za zgornji primer verjetno nekaj v slogu:
a = [1; 3; 6; 7; 10; 12];
b = [2; 3; 8; 17; 25; 44];

[p, q] = meshgrid(a, b);
r = [p(:) q(:)];

for i = r(:,1)
    for j = r(:,2)
        X = i + j;
    end
end

TEDY ::

če sta a in b vektroja ne potrebuje for zank.
a = [1, 3, 6, 7, 10, 12]
b = [2, 3, 8, 17, 25, 44]

X = a' + b
X = a' .* b
X = bsxfun(@plus, a', b) % za poljubne operatorje
...

Zgodovina sprememb…

  • spremenil: TEDY ()

kitajc ::

Hvla vsem, deluje pa najbolje tole (brez ' pred a-ji v 5. in 6. vrstici). Matlab uporabljam od petka in dejansko rabim zelo podrobna navodila.

Matrika deluje le na enako število različnih vrednosti a in b. Sam imam pri a 15 vrednosti, pri b pa 50 vrednosti. Ok, to bom le pobrisal zadnje stolpce.

Kako pa sedaj v kodo vstavim enačbo, ki ima notri vrednosti, ki jih pobira od prej.

Npr.: X = a + b + T + 20

T sem določil že v prejšnjih enačbah in je le ENO število. Npr. T = 3.

T=3
a = [1, 3, 6, 7]
b = [2, 3, 8, 17]

X = a + b + T
X = a .* b + T
X = bsxfun(@plus, a', b) % za poljubne operatorje


Vendar dobim tole:

T =

3


a =

1 3 6 7


b =

2 3 8 17


X =

6 9 17 27


X =

5 12 51 122


X =

3 4 9 18
5 6 11 20
8 9 14 23
9 10 15 24



BTW - moja končna enačba pa je taka (vstavljam različne vrednosti Pr_0 in Tr_0, ostalo dobim od prej)
syms Vr_0
eqn = Pr_0 * Vr_0/Tr_0 == 1 + B_0/Vr_0 + C_0/Vr_0^2 + D_0/Vr_0^5 + c4_0/(Tr_0^3 * Vr_0^2) * (BETA_0 + (GAMA_0/Vr_0^2)* exp(-GAMA_0/Vr_0^2)) ;
Vr_0 = double(solve(eqn))

TEDY ::

V Octave moj primer deluje...

Vzami primer of bluefish in spremeni
X = i + j;
v
X = i + j + T + 20;


vendar se v matlabu probaj čimbolj izogibati zankam (for, while)...
zato lahko uporabiš
X = p(:) + q(:) + T + 20

namesto dveh for zank

kitajc ::

Za tole kodo:

T = 19
a = [1; 3; 6; 7];
b = [2; 3; 8; 17];
[p, q] = meshgrid(Pr_0, Tr_0);
r = [p(:) q(:)];
X = p(:) + q(:) + T + 20
X = bsxfun(@plus, a', b)

Dobim:

T =

19


X =

42
43
48
57
44
45
50
59
47
48
53
62
48
49
54
63


X =

3 5 8 9
4 6 9 10
9 11 14 15
18 20 23 24


Vrednosti X-a se spreminjajo, če spreminjam T, matrika pa ostane vsekozi enaka. Ne zastopim zadnjih dveh vrstic. Kako je tu sploh vstavljena enačba oz. kako bom sam vstavil svojo?

Zgodovina sprememb…

  • spremenilo: kitajc ()

kitajc ::

Recimo kasneje bom imel tole enačbo. Odebeljene vrednosti so izračunane že od prej. Podčrtane vstavljam. Končni rezultat bo Vr.

Pr_0 * Vr_0/Tr_0 == 1 + B_0/Vr_0 + C_0/Vr_0^2 + D_0/Vr_0^5 + c4_0/(Tr_0^3 * Vr_0^2) * (BETA_0 + (GAMA_0/Vr_0^2)* exp(-GAMA_0/Vr_0^2))

TEDY ::

zadnja vrstica je odveč, več ti pa žal ne znam pomagati

chort ::

Kot je rekel TEDY, zadnja vrstica je odveč (oz ni odveč, le ne prišteva še T-ja in fiksnega števila 20) Sicer je pa prvi X ki ga dobiš pravilen.

Končno enačbo premeči tako, da boš imel na levi samo še Vr_0 (Vr_0 == klobasa na desni * ( Tr_0/Pr_0) ), preveri mojo logiko ker sem iz vaje.

Poskusi kaj dobiš če v desni strani premetane enačbe tvoje podčrtane parametre zamenjaš z ustreznimi klici na vektorje Pr_0(:), Tr_(:). Če dobiš preveč rezultatov (ker izračunalo kombinacijo 3 spremenljivk - Tr_0 boš imel na desni strani dvakrat), se vrni na rešitev ki ti jo je podal Bluefish in tisti njegov X = i + j zamenjaj s svojo formulo.

kitajc ::

To sem rešil takole. Koda se poračuna in dobim matriko.

Tr_0=[1.34, 1.515789]
Pr_0=[0.0220319, 1.345, 1.5]

for j =1:2
for i = 1:3
syms Vr_0
solve(-Pr_0(i)*Vr_0/Tr_0(j) + 1 + B_0/Vr_0 + C_0/Vr_0^2 + D_0/Vr_0^5 + c4_0/(Tr_0(j)^3*Vr_0^2)*(BETA_0 + GAMA_0/Vr_0^2)*exp(-GAMA_0/Vr_0^2), Vr_0)
end
end
Vr_0 = bsxfun(@plus, Tr_0', Pr_0)


Sedaj pa bi moral vse te rezultate za Vr_0 vstaviti v končno enačbo za Z0_0 in dobiti tabelo rezultatov. Poskusil sem z še eno for zanko, a ne gre. Če kodo dodam na koncu, pa tudi ne. Drugi del kode je podoben prvemu in izgleda takole:
%Izračun Z0_0
for k =1:6
syms Z0_0
solve(-Z0_0 + 1 + B_0/Vr_0(k) + C_0/Vr_0(k)^2 + D_0/Vr_0(k)^5 + c4_0/(Tr_0^3*Vr_0(k)^2)*(BETA_0 + (GAMA_0/Vr_0(k)^2))*exp(-GAMA_0/Vr_0(k)^2), Z0_0)
end
Z0_0 = bsxfun(@plus, Tr_0', Pr_0)


Kaj lahko storim? Zakaj mi ne deluje, saj sta kodi lepo zaporedno razporejeni?

kitajc ::

Tole je, kar mi je uspelo do sedaj. V to enačbo sem vstavljal Tr_0 in Pr_0.
Tr_0=[1, 1.5]
Pr_0=[0.5, 1.3, 2]

for j =1:2
for i = 1:3
syms Vr_0
solve(-Pr_0(i)*Vr_0/Tr_0(j) + 1 + B_0/Vr_0 + C_0/Vr_0^2 + D_0/Vr_0^5 + c4_0/(Tr_0(j)^3*Vr_0^2)*(BETA_0 + GAMA_0/Vr_0^2)*exp(-GAMA_0/Vr_0^2), Vr_0)
end
end
Vr_0 = bsxfun(@plus, Tr_0', Pr_0)


Končni rezultat je matrika:
Vr_0 =

1.362031900000000 2.685000000000000 2.840000000000000
1.537820900000000 2.860789000000000 3.015789000000000


Sedaj pa moram v naslednjo enačbo vstaviti Vr_0 (rezultati matrike) in spet Tr_0 iz prvega dela. Dobiti moram neko matriko, da bo lepo pregledno.
%Izračun Z0_0
for j = Vr_0
for i = 1:3
syms Z0_0
solve(-Z0_0 + 1 + B_0/Vr_0 + C_0/Vr_0^2 + D_0/Vr_0^5 + c4_0/(Tr_0^3*Vr_0^2)*(BETA_0 + (GAMA_0/Vr_0^2))*exp(-GAMA_0/Vr_0^2), Z0_0)
end
end
Z0_0 = bsxfun(@plus, Tr_0', Pr_0)

Na tak način, kot zgoraj... Kako sploh pokličem matriko?

TEDY ::

v tvojem primeru, ko kličeš solve() nikamor ne shraniš rezultata.

%clear all
Tr_0=[1.34, 1.515789]
Pr_0=[0.0220319, 1.345, 1.5]
syms Vr_0 Z0_0

%B_0 = 1; C_0 = 1; D_0 = 1; c4_0 = 1; BETA_0 = 1; GAMA_0 = 1;

Z0 = [];
for j =1:length(Tr_0)
    for i = 1:length(Pr_0)
        Vr = double(solve(-Pr_0(i)*Vr_0/Tr_0(j) + 1 + B_0/Vr_0 + C_0/Vr_0^2 + D_0/Vr_0^5 + c4_0/(Tr_0(j)^3*Vr_0^2)*(BETA_0 + GAMA_0/Vr_0^2)*exp(-GAMA_0/Vr_0^2), Vr_0));
        
        %Izračun Z0_0
        Z0(end+1) = double(solve(-Z0_0 + 1 + B_0/Vr + C_0/Vr^2 + D_0/Vr^5 + c4_0/(Tr_0(j)^3*Vr^2)*(BETA_0 + (GAMA_0/Vr^2))*exp(-GAMA_0/Vr^2), Z0_0));
    end
end

disp(Z0)

kitajc ::

Res? Ker jih v Workspace vseskozi vidim. Tudi izračune za Vr_0, ki sem jih omenil zgoraj. V Workspace-u 2x kliknem na matriko in odpre mi tabelo rezultatov (tako kot v Excelu).

Kako pa bi lahko tole transcendentno enačbo rešil drugače?

Oz. kako rešim drugi del, kjer vstavljam - kličem matriko Vr_0 in vrednosti Tr_0 (tiste od 0.3 pa do 4.00)?

TEDY ::

vidiš rezultate (žal napačne), ker ti jih izračuna tole Vr_0 = bsxfun(@plus, Tr_0', Pr_0) takoj za for zankami.
Zgoraj sem ti dal primer, ko se rezultat shrani v Vr, ta pa se nato uporabi za izračun Z0

kitajc ::

Tale koda deluje odlično. Teh funkcij ne poznam in ne vem, kaj počnem. Sem se pa malo zafrknil, ker sem prej računal te kontante B, C in D le za eno spremenljivko. Sedaj bi moral v kodo vstaviti enačbe vseh treh, ki pa nabira vrednosti le Tr_0.

Verjetno se to da narediti s temi zankami? Vse tiste b1, b2 itd pa so standardi. Tele 3 pa je treba izračunati vsakič posebej.

B = double(solve(-B_0 + b1_0 - b2_0/Tr_0(j) - b3_0/Tr_0(j)^2 - b4_0/Tr_0(j)^3))
C = double(solve(-C_0 + c1_0 - c2_0/Tr_0(j) + c3_0/Tr_0(j)^3))
D = double(solve(-D_0 + d1_0 + d2_0/Tr_0(j)))


Ali sem na pravi poti s spodnjo kodo?

Z0 = []; 
    for j = 1:length(Tr_0)
        B = double(solve(-B_0 + b1_0 - b2_0/Tr_0(j) - b3_0/Tr_0(j)^2 - b4_0/Tr_0(j)^3))
        C = double(solve(-C_0 + c1_0 - c2_0/Tr_0(j) + c3_0/Tr_0(j)^3))
        D = double(solve(-D_0 + d1_0 + d2_0/Tr_0(j)))
    for i = 1:length(Pr_0)
            Vr = double(solve(-Pr_0(i)*Vr_0/Tr_0(j) + 1 + B_0/Vr_0 + C_0/Vr_0^2 + D_0/Vr_0^5 + c4_0/(Tr_0(j)^3*Vr_0^2)*(BETA_0 + GAMA_0/Vr_0^2)*exp(-GAMA_0/Vr_0^2), Vr_0));
         
        %Izračun Z0_0
        Z0(end+1) = double(solve(-Z0_0 + 1 + B_0/Vr + C_0/Vr^2 + D_0/Vr^5 + c4_0/(Tr_0(j)^3*Vr^2)*(BETA_0 + (GAMA_0/Vr^2))*exp(-GAMA_0/Vr^2), Z0_0));
    end
end
 
disp(Z0)


Kaj pomeni (end+1)?

TEDY ::

Si na pravi poti, vendar na novo izračunanih B,C,D ne upoštevaš v spodnjih dveh enačbah. Če so vse spremenljivke znane lahko poenostaviš:
B = -B_0 + b1_0 - b2_0/Tr_0(j) - b3_0/Tr_0(j)^2 - b4_0/Tr_0(j)^3;
C = -C_0 + c1_0 - c2_0/Tr_0(j) + c3_0/Tr_0(j)^3;
D = -D_0 + d1_0 + d2_0/Tr_0(j);

(end) ti vrne index zadnjega elementa v vektorju Z0, ker pa nočeš, da ti razultat prepiše zadnji element ga zapiše na +1 mesto naprej.

Zgodovina sprememb…

  • spremenil: TEDY ()

kitajc ::

%Konstante za osnovni fluid
b1_0=0.1181193; b2_0=0.265728; b3_0=0.154790; b4_0=0.03023; c1_0=0.0236744; c2_0=0.0186984; c3_0=0; c4_0=0.042724; d1_0=0.155488*10^-4; 
d2_0=0.623689*10^-4; BETA_0=0.65392; GAMA_0=0.060167;

%clear all
Tr_0=[1.34, 1.515789]
Pr_0=[0.0220319, 1.345, 1.5]
syms B_0 C_0 D_0 Vr_0 Z0_0
 
Z0 = []; 
    for j = 1:length(Tr_0)
        B = (-B_0 + b1_0 - b2_0/Tr_0(j) - b3_0/Tr_0(j)^2 - b4_0/Tr_0(j)^3);
        C = (-C_0 + c1_0 - c2_0/Tr_0(j) + c3_0/Tr_0(j)^3);
        D = (-D_0 + d1_0 + d2_0/Tr_0(j));
    for i = 1:length(Pr_0)
            Vr = double(solve(-Pr_0(i)*Vr_0/Tr_0(j) + 1 + B_0/Vr_0 + C_0/Vr_0^2 + D_0/Vr_0^5 + c4_0/(Tr_0(j)^3*Vr_0^2)*(BETA_0 + GAMA_0/Vr_0^2)*exp(-GAMA_0/Vr_0^2), Vr_0));
          
        %Izračun Z0_0
        Z0(end+1) = double(solve(-Z0_0 + 1 + B_0/Vr + C_0/Vr^2 + D_0/Vr^5 + c4_0/(Tr_0(j)^3*Vr^2)*(BETA_0 + (GAMA_0/Vr^2))*exp(-GAMA_0/Vr^2), Z0_0));
    end
end
 disp(Z0)


Error:
Error using  ^ 
Inputs must be a scalar and a square matrix.
To compute elementwise POWER, use POWER (.^) instead.

Error in dfgherthe (line 51)
        Z0(end+1) = double(solve(-Z0_0 + 1 + B_0/Vr + C_0/Vr^2 + D_0/Vr^5
        + c4_0/(Tr_0(j)^3*Vr^2)*(BETA_0 +
        (GAMA_0/Vr^2))*exp(-GAMA_0/Vr^2), Z0_0));

PRINCIP NAJ BI BIL TAK:
1. Konstante (b1, b2...gama, beta) in Tr (od 0.3 pa do 4) vstavim v enačbe za B, C in D.
2. B, C, D + Tr (od 0.3 pa do 4) + Pr (od 0.05 do 10) vstavljam za enačbo Vr in dobim matriko rezultatov za Vr. Matrika, kjer so na y osi vstavljene vrednosti Tr, na x pa Pr.
3. Vse vstavim v zadnjo enačbo za Z0, kjer dobim matriko rezultatov Z0. Matrika, kjer so na y osi vstavljene vrednosti Tr, na x pa Pr. (to je edina matrika, ki jo dejansko rabim videti)

Je pa tole postalo zame malo prekomplicirano. Upam, da bo to na koncu vse zmešano.

TEDY ::

%Konstante za osnovni fluid
b1_0=0.1181193; b2_0=0.265728; b3_0=0.154790; b4_0=0.03023; c1_0=0.0236744; c2_0=0.0186984; c3_0=0; c4_0=0.042724; d1_0=0.155488*10^-4; 
d2_0=0.623689*10^-4; BETA_0=0.65392; GAMA_0=0.060167;
 
%clear all
Tr_0=[1.34, 1.515789]
Pr_0=[0.0220319, 1.345, 1.5]
syms Vr_0 Z0_0
  
Z0 = []; 
    for j = 1:length(Tr_0)
        B_0 = (b1_0 - b2_0/Tr_0(j) - b3_0/Tr_0(j)^2 - b4_0/Tr_0(j)^3);
        C_0 = (c1_0 - c2_0/Tr_0(j) + c3_0/Tr_0(j)^3);
        D_0 = (d1_0 + d2_0/Tr_0(j));
    for i = 1:length(Pr_0)
            Vr = double(solve(-Pr_0(i)*Vr_0/Tr_0(j) + 1 + B_0/Vr_0 + C_0/Vr_0^2 + D_0/Vr_0^5 + c4_0/(Tr_0(j)^3*Vr_0^2)*(BETA_0 + GAMA_0/Vr_0^2)*exp(-GAMA_0/Vr_0^2), Vr_0));
           
        %Izračun Z0_0
        Z0(end+1) = double(solve(-Z0_0 + 1 + B_0/Vr + C_0/Vr^2 + D_0/Vr^5 + c4_0/(Tr_0(j)^3*Vr^2)*(BETA_0 + (GAMA_0/Vr^2))*exp(-GAMA_0/Vr^2), Z0_0));
    end
end
 disp(Z0)

Probaj tole, jaz trenutno ne morem sprobati.

kitajc ::

Sem poskusil in deluje, vendar le za določene vrednosti Tr in Pr. Verjetno je fora v tisti transcendenti enačbi za Vr. V nekaterih delih tabele dobim pravilne rezultate (to velja za 85-90% tabele). Ko pa se pomikam proti manjšim vrednostim Tr (< 0,5) in Pr (< 0,4), pa so le še nekateri rezultati pravilni, na koncu pa so vsi povsem napačni.

Bolj pa me zanima naslednje. V nadaljevanju bom moral uporabiti to tabelo Z-jev za izračun različnih veličin. Tam bodo nastopali le Z-ji ter Tr in Pr (isti kot prej).

Npr.: Tole je enačba za izračun fugativnosti. Vsakič bom moral poračunati E (vseskozi enaka enačba) in nato vstaviti v veličino (5 jih imam, ena ima še nekaj zraven ostale so na isti način). To so potem končni rezultati. Vstavljati rabim dobljene vrednosti Vr in Z.

fugativnost = ln (f/p) = Z - 1 - ln(Z) + B/(2*Vr^2) + D/(5*Vr^5) + E

E = c4/(2*Tr^3*GAMA) (BETA + 1 - (BETA + 1 + GAMA/Vr^2)*exp(-GAMA/Vr^2))


Najbrž bi moral potem odpreti nov .m dokument, poklicati vrednosti Z-jev ali pa Vr-jev in potem napisati kodo za izračun npr. fugativnosti. Če ne bom znal rešiti tistega prvega problema, potem bi to lahko naredil za del, ki deluje avtomatično, ostalo pa spet poračunal kako drugače. Ali bom lahko iz te kode klical tudi Vr-je?

Druače je pa tole super koda!

Zgodovina sprememb…

  • spremenilo: kitajc ()

kitajc ::

Dobim prvih 47 rezultatov, potem pa sledi:

Warning: Cannot find explicit solution.
> In solve at 319
In izracun_z0 at 17
Error using ^
Inputs must be a scalar and a square matrix.
To compute elementwise POWER, use POWER (.^) instead.

Error in izracun_z0 (line 20)
Z0(end+1) = double(solve(-Z0_0 + 1 + B_0/Vr + C_0/Vr^2 + D_0/Vr^5
+ c4_0/(Tr_0(j)^3*Vr^2)*(BETA_0 +
(GAMA_0/Vr^2))*exp(-GAMA_0/Vr^2), Z0_0));


Moram pogledati še, če v teoriji piše kaj od teh majhnih vrednosti...

Ko sem poskusil vstavljati le nekatere manjše vrednosti Tr in Pr, pa sem VEDNO dobil rezultate, a žal napačne.

Ne, ne dobim vedno!

Zgodovina sprememb…

  • spremenilo: kitajc ()

kitajc ::

Kako nastavim spodnjo kodo? Izračunane B, C in D vstavljam naprej v enačbi za Vr in Z0, kar je bilo prej v redu. Sedaj sem dodal še da se rezultati vstavijo v E in nato v F0, a mi jih ne bere.

%REDUCIRANA FUGATIVNOST ln(f/p)
syms Vr_0 Z0_0 E F0
F0 = []; 
    for j = 1:length(Tr_0)
        B_0 = (b1_0 - b2_0/Tr_0(j) - b3_0/Tr_0(j)^2 - b4_0/Tr_0(j)^3);
        C_0 = (c1_0 - c2_0/Tr_0(j) + c3_0/Tr_0(j)^3);
        D_0 = (d1_0 + d2_0/Tr_0(j));
    for i = 1:length(Pr_0)
            Vr = double(solve(-Pr_0(i)*Vr_0/Tr_0(j) + 1 + B_0/Vr_0 + C_0/Vr_0^2 + D_0/Vr_0^5 + c4_0/(Tr_0(j)^3*Vr_0^2)*(BETA_0 + GAMA_0/Vr_0^2)*exp(-GAMA_0/Vr_0^2), Vr_0));
            Z0 = double(solve(-Z0_0 + 1 + B_0/Vr + C_0/Vr^2 + D_0/Vr^5 + c4_0/(Tr_0(j)^3*Vr^2)*(BETA_0 + (GAMA_0/Vr^2))*exp(-GAMA_0/Vr^2), Z0_0));
            E = double(solve(-E + c4_0/(2*Tr_0^3) * (BETA + 1 - (BETA + 1 + GAMA/Vr^2)* exp(-GAMA/Vr^2))), E);
        %Izračun F0
        F0(end+1) = double(solve(-Z0_0 + 1 + B_0/Vr + C_0/Vr^2 + D_0/Vr^5 + c4_0/(Tr_0(j)^3*Vr^2)*(BETA_0 + (GAMA_0/Vr^2))*exp(-GAMA_0/Vr^2), F0));
    end 
end
disp(F0)

TEDY ::

v E bi moral vstaviti Tr_0(j) namesto Tr_0 ?
v F0 bi moral vstaviti Z0 namesto Z0_0 ?

kitajc ::

Res je, pa tudi tiste GAME in BETE niso bile ok.

Malo sem popravil, sedaj pa dobim nek čuden error (line 40 je v tem zapisu vrstica 18):
Error using sym/double
Too many input arguments.

Error in fugativnost (line 40)
            E = double(solve(-E + c4_0/(2*Tr_0(j)^3) * (BETA_0 + 1 -
            (BETA_0 + 1 + GAMA_0/Vr^2)* exp(-GAMA_0/Vr^2))), E);


Koda je taka:
b1_0=0.1181193; b2_0=0.265728; b3_0=0.154790; b4_0=0.03023; c1_0=0.0236744; c2_0=0.0186984; c3_0=0; c4_0=0.042724; d1_0=0.155488*10^-4; 
d2_0=0.623689*10^-4; BETA_0=0.65392; GAMA_0=0.060167;

%clear all
Tr_0=[0.85, 0.90, 0.93, 0.95]
Pr_0=[0.800, 1.000, 1.200, 1.500]

%REDUCIRANA FUGATIVNOST ln(f/p)
syms Vr_0 Z0_0 E F0
F0 = []; 
    for j = 1:length(Tr_0)
        B_0 = (b1_0 - b2_0/Tr_0(j) - b3_0/Tr_0(j)^2 - b4_0/Tr_0(j)^3);
        C_0 = (c1_0 - c2_0/Tr_0(j) + c3_0/Tr_0(j)^3);
        D_0 = (d1_0 + d2_0/Tr_0(j));
    for i = 1:length(Pr_0)
            Vr = double(solve(-Pr_0(i)*Vr_0/Tr_0(j) + 1 + B_0/Vr_0 + C_0/Vr_0^2 + D_0/Vr_0^5 + c4_0/(Tr_0(j)^3*Vr_0^2)*(BETA_0 + GAMA_0/Vr_0^2)*exp(-GAMA_0/Vr_0^2), Vr_0));
            Z0 = double(solve(-Z0_0 + 1 + B_0/Vr + C_0/Vr^2 + D_0/Vr^5 + c4_0/(Tr_0(j)^3*Vr^2)*(BETA_0 + (GAMA_0/Vr^2))*exp(-GAMA_0/Vr^2), Z0_0));
            E = double(solve(-E + c4_0/(2*Tr_0(j)^3) * (BETA_0 + 1 - (BETA_0 + 1 + GAMA_0/Vr^2)* exp(-GAMA_0/Vr^2))), E);
        %Izračun F0
        F0(end+1) = double(solve(-Z0 + 1 + B_0/Vr + C_0/Vr^2 + D_0/Vr^5 + c4_0/(Tr_0(j)^3*Vr^2)*(BETA_0 + (GAMA_0/Vr^2))*exp(-GAMA_0/Vr^2), F0));
    end 
end
disp(F0)

TEDY ::

čist tazadnji "E" si en oklepaj preveč ven dal...

kitajc ::

Vidim, na koncu morata biti dva. Tole sem dal pa na isti način kot prej, a mi javi:

Error using solve>getEqns (line 422)
The input contains an empty equation or variable.

Error in solve (line 227)
[eqns,vars,options] = getEqns(varargin{:});

Error in fugativnost (line 42)
            F0 (end+1) = double(solve(-F0 + Z0 + 1 - log(Z0) + B_0/Vr +
            C_0/(2*Vr^2) + D_0/(5*Vr^5) + E, F0));


%REDUCIRANA FUGATIVNOST ln(f/p)
syms Vr_0 Z0_0 E F0
F0 = []; 
    for j = 1:length(Tr_0)
        B_0 = (b1_0 - b2_0/Tr_0(j) - b3_0/Tr_0(j)^2 - b4_0/Tr_0(j)^3);
        C_0 = (c1_0 - c2_0/Tr_0(j) + c3_0/Tr_0(j)^3);
        D_0 = (d1_0 + d2_0/Tr_0(j));
    for i = 1:length(Pr_0)
            Vr = double(solve(-Pr_0(i)*Vr_0/Tr_0(j) + 1 + B_0/Vr_0 + C_0/Vr_0^2 + D_0/Vr_0^5 + c4_0/(Tr_0(j)^3*Vr_0^2)*(BETA_0 + GAMA_0/Vr_0^2)*exp(-GAMA_0/Vr_0^2), Vr_0));
            Z0 = double(solve(-Z0_0 + 1 + B_0/Vr + C_0/Vr^2 + D_0/Vr^5 + c4_0/(Tr_0(j)^3*Vr^2)*(BETA_0 + (GAMA_0/Vr^2))*exp(-GAMA_0/Vr^2), Z0_0));
            E = double(solve(-E + c4_0/(2*Tr_0(j)^3) * (BETA_0 + 1 - (BETA_0 + 1 + GAMA_0/Vr^2)* exp(-GAMA_0/Vr^2)), E));
            %Izračun F0
            F0 (end+1) = double(solve(-F0 + Z0 + 1 - log(Z0) + B_0/Vr + C_0/(2*Vr^2) + D_0/(5*Vr^5) + E, F0));
    end
end
disp(F0)


Sedaj se s temle ubadam že nekaj časa. A je kaj zaradi tega (end+1)?

Je pa v Matlabu funkcija ln enaka kar log? Tam bi sicer moralo biti ln(Z0).

Zgodovina sprememb…

  • spremenilo: kitajc ()

kitajc ::

Sem nekako rešil:

%REDUCIRANA FUGATIVNOST ln(f/p)
syms Vr_0 Z0_0 E0 F
F0 = []; 
    for j = 1:length(Tr_0)
        B_0 = (b1_0 - b2_0/Tr_0(j) - b3_0/Tr_0(j)^2 - b4_0/Tr_0(j)^3);
        C_0 = (c1_0 - c2_0/Tr_0(j) + c3_0/Tr_0(j)^3);
        D_0 = (d1_0 + d2_0/Tr_0(j));
    for i = 1:length(Pr_0)
            Vr = double(solve(-Pr_0(i)*Vr_0/Tr_0(j) + 1 + B_0/Vr_0 + C_0/Vr_0^2 + D_0/Vr_0^5 + c4_0/(Tr_0(j)^3*Vr_0^2)*(BETA_0 + GAMA_0/Vr_0^2)*exp(-GAMA_0/Vr_0^2), Vr_0));
            Z0 = double(solve(-Z0_0 + 1 + B_0/Vr + C_0/Vr^2 + D_0/Vr^5 + c4_0/(Tr_0(j)^3*Vr^2)*(BETA_0 + (GAMA_0/Vr^2))*exp(-GAMA_0/Vr^2), Z0_0));
            E = double(solve(-E0 + c4_0/(2*Tr_0(j)^3) * (BETA_0 + 1 - (BETA_0 + 1 + GAMA_0/Vr^2)* exp(-GAMA_0/Vr^2)), E0));
            %Izračun F0
            F0 (end+1) = double(solve(-F + Z0 + 1 - log(Z0) + B_0/Vr + C_0/(2*Vr^2) + D_0/(5*Vr^5) + E, F));
    end
end
disp(F0)

bluefish ::

Ja, v Matlabu je log uporabljen za naravni logaritem, log10 pa za desetiškega.


Vredno ogleda ...

TemaSporočilaOglediZadnje sporočilo
TemaSporočilaOglediZadnje sporočilo
»

Thinkpad t420 WIN10

Oddelek: Pomoč in nasveti
403210 (2259) Klemzz
»

Kako poračunati enačbo?

Oddelek: Loža
101544 (1063) kitajc
»

matematka

Oddelek: Šola
233030 (2009) lebdim
»

počasno delovanje kompa Asus N76VZ

Oddelek: Pomoč in nasveti
81090 (918) WIngs
»

BSOD me spravlja v obup

Oddelek: Strojna oprema
251875 (1347) rsroki

Več podobnih tem