Winkel iterieren mit Scilab

Neue Frage »

stevo-frbg Auf diesen Beitrag antworten »
Winkel iterieren mit Scilab
Meine Frage:
Hallo liebes Matheboard,

ich schreibe gerade eine Arbeit über Windenergie, und da soll ein Winkel bestimmt werden. Problem ist, die gegebene Gleichung ist nicht einfach auflösbar, sondern muss iterativ gelöst werden.

Wie mache ich das am besten mit Scilab? Ich habe schon ein bisschen rumgebastelt, klappt aber nicht so richtig.

Die Gleichung lautet:

t * cA(alpha - alpha_p) - (8*pi*r)/z *sin(alpha) * tan(alpha_1-alpha) = 0

t, cA-Kurve, alpha_p, r, z und alpha_1 sind bekannt, bestimmt werden soll also alpha.

Wie löse ich sowas denn am elegantesten? Theoretisch hat man das alles ja mal im Grundstudium gelernt, aber irgendwie ist es alles so lange her... :-(

Danke schon jetzt für eure Mühe!

Gruß, Stevo

Meine Ideen:
Ich habe versucht, das ganze per while- schleife durchzuvariieren, aber irgendwie rechnet der PC sich da zu Tode.

Das bisherige Programm sieht, gekürzt aufs wesentliche, so aus:
----------------------------------------------------------------

function [alpha,i] = iteration(t,z,r,R,alpha_1,alpha_pitch)

i = 1
f=1

while(abs(f(i))>=0.1)


alpha_A = alpha - alpha_pitch


exec('C:\Users\steffen\Documents\Studium\Thesis\Blattelementmethode\ca_cw_f
unc.sce', -1)
[ca,cw] = ck220cacw(alpha_A) //Getestet, funktioniert

y = x/sin(2/3*alpha_1)

f(i) = t*ca - 8*%pi*r/z *sin(alpha)*tan(alpha_1-alpha)

if f(i)>0 then alpha = alpha*0.95
//Prüfung: alpha zu klein oder zu groß?
elseif f(i)<0 then alpha = alpha*1.05
end

i=i+1

end

endfunction
Steffen Bühler Auf diesen Beitrag antworten »
RE: Winkel iterieren mit Scilab
Hier bietet sich die binäre Suche an. Du hast ein erstes Intervall für alpha. Da berechnest Du für die beiden Grenzen die Funktionswerte. Ist einer positiv und der andere negativ, muß dazwischen die Nullstelle liegen. Nun berechne den Funktionswert der Intervallmitte. Je nach dessen Polarität kannst Du nun in der oberen oder unteren Intervallhälfte weitermachen, solange, bis Du nah genug an die Nullstelle herangekommen bist.

Viele Grüße
Steffen
stevo-frbg Auf diesen Beitrag antworten »

Hallo Steffen,

vielen Dank für deinen Tipp.

Gruß,
Steffen
stevo-frbg Auf diesen Beitrag antworten »

Hallo,

irgendwie klappt das nicht.

Die Iteration lautet jetzt so:
if f > 0 then o = alpha //Prüfung: alpha zu klein oder zu groß?
elseif f < 0 then u = alpha //Prüfung nach der BINÄREN ITERATION
elseif abs(f) <= 0.002 then break;
end

Problem scheint zu sein, dass oft schon Obergrenze und Untergrenze ein positives f-ergeben. Gibt es vielleicht noch eine andere Methode?

Ich habe eben noch versucht, einfach einen Vektor für alpha zu definieren und den einzusetzen, und dann per ab(min(f)) den Wert zu finden, der am nächsten bei 0 liegt. Aber selbst das funktioniert nicht wirklich zufriedenstellend. :-(
Steffen Bühler Auf diesen Beitrag antworten »

Diese Methode setzt aber monotones Steigen voraus, f kann ja auch fallen. Weißt Du einigermaßen, wie f aussieht? Du solltest schon zwei alphas kennen, zwischen denen die Nullstelle liegt. Ganz ohne geht's nicht.

Wenn es eine Nullstelle gibt, muß es auch zwei Werte links und rechts davon geben, die unterschiedliches Vorzeichen haben. Das ist das Ursprungsintervall, bei Dir begrenzt durch u und o. In Deinen zwei Zeilen sehe ich jetzt nicht, wie dieses Intervall halbiert wird, aber es sollte ungefähr so laufen:

code:
1:
2:
3:
4:
5:
6:
7:
8:
9:
10:
11:
12:
13:
14:
15:
16:
17:
18:
19:
20:
21:
22:
23:
u = -10               'zum Beispiel
o = 10                'zum Beispiel
fmitte = 1            'damit es erstmal losläuft

while (abs(fmitte) >= 0.002)

  mitte = (u+o)/2

  alpha = u
  fu = t*ca - 8*%pi*r/z *sin(alpha)*tan(alpha_1-alpha)
  alpha = o
  fo = t*ca - 8*%pi*r/z *sin(alpha)*tan(alpha_1-alpha)
  alpha = mitte
  fmitte = t*ca - 8*%pi*r/z *sin(alpha)*tan(alpha_1-alpha)

  if sgn(fu) = sgn(fmitte)      'nullstelle liegt nicht in der unteren Hälfte
    u = mitte                   'untere Grenze rauf

  else if sgn(fo) = sgn(fmitte) 'nullstelle liegt nicht in der oberen Hälfte
    o = mitte                   'obere Grenze runter

wend

Da könnten Syntaxfehler drin sein, mein Basic ist ein wenig eingerostet. Aber ich glaube, es ist klar, was ich meine.

Viele Grüße
Steffen
stevo-frbg Auf diesen Beitrag antworten »

Irgendwie klappt das einfach nicht, die Funktion der Winkel sieht am Schluss so aus:

(siehe Anhang)

Die maximale Iterationszahl von zB. 500 wird auch fast jedes mal erreicht - das Ergebnis scheint also extrem daneben zu liegen!

Ich habe nun versucht, das ganze per fsolve-Befehl, welcher ja per Newton-Verfahren löst, klarzukommen - ein Freund meinte, das Binärverfahren wäre nicht so praktisch.

Ich kriege es aber garnicht hin - Aufruf lautet bisher:

------------------------

exec('C:\Users\steffen\Documents\Studium\Thesis\Blattelementmethode\fx63-137.sce', -1) //Ruft das Programm auf, welches den dazu passenden ca-Wert ausgibt - müsste eigentlich stimmen!

R = 3.6
dr = 0.1
rmin = 0.1
r = rmin:dr:R


z=3
t=0.22

lambda = 7
alfa_pitch = 0
deff('[f]=gleichung(x)','f = t*ca - 8*%pi*r_loc/z *sin(x)*tan(alfa_1-x)')

for i=1:length(r)

alfa_1 = atan(R/lambda*r(i))
alfa_A = alfa_1-alfa_pitch
[ca,cw] = cacw(alfa_A)
r_loc = r(i)
alfa(i)=fsolve(alfa_1,gleichung)

end

plot(r,alfa)


Leider ist das Ergebnis garnicht plausibel (Winkel wächst mit steigendem r an, anstatt, wie es sein sollte, von etwa 90° ab abzufallen - da stimmt die ursprüngliche Lösung tendenziell schon!)

Viele Grüße
 
 
Neue Frage »
Antworten »



Verwandte Themen

Die Beliebtesten »
Die Größten »
Die Neuesten »