Temperaturumrechnung PT1000

Für sonstige Unterhaltungen, welche nicht direkt mit Lazarus zu tun haben
Antworten
siro
Beiträge: 730
Registriert: Di 23. Aug 2016, 14:25
OS, Lazarus, FPC: Windows 11
CPU-Target: 64Bit
Wohnort: Berlin

Temperaturumrechnung PT1000

Beitrag von siro »

Hallo zusammen,
eigentlich wollte ich "nur mal eben" die Formel umstellen für einen Temperatursensor PT100.
nachdem ich 2 Bekannte mit "wesentlich" besseren Mathekenntnissen um Hilfe gebeten habe und beide meinten,
das es ihnen zu viel Arbeit ist und zudem es wohl bis zu 4 Lösungen gibt, bin ich "etwas" verwundert.
Okay, das liegt wohl daran, dass meine Mathekenntnisse unterstes Niveau sind.
Aber ich hätte nicht gedacht, dass es so kompliziert sein kann...

Ich möchte auf diesem Wege hier bitte niemanden mit komplizierten Aufgaben versorgen,
aber eventuell, hat das ja schon mal jemand gemacht oder hat einen Link wo die Formel schon mal umgestellt wurde.

gegeben:
Widerstand(R) = R0 *(1 + A*t + B*t*t + C*(t-100)*t*t*t); // diese Formel möchte ich umgestellt haben nach "t"

Es handelt sich wohl um eine "Quartische Gleichung" Polynom 4ten Grades
https://de.wikipedia.org/wiki/Quartische_Gleichung

Code: Alles auswählen

 
Mein momentaner Pascal Code sieht so aus:
Für die positiven Temperaturen "Function CalcT" läuft schon.
 
var R0 : Double =  1000;       // Widerstand in Ohm bei 0 Grad Celsius Wert für PT1000
var: Double =  3.9083E-3;  // Koeffizient 1
var: Double = -5.775E-7;   // Koeffizient 2
var: Double = -4.183E-12;  // Koeffizient 3
 
//------------------------------------------------------------------------------
// Ermittelt die Temperatur bei dem übergebenen Widerstandswert R
// !!!! Formel momentan nur für positive Temperaturn gueltig !!!!
// Bei negativen Temperaturen, also Widerstandswerte unter R0,
// wird die Berechnung anscheinend wesentlich komplexer....
 
 
Function CalcT(R:Double):Double;
var x,y:Double;
begin
  x := A / (B*2);                  // Zwischenvariable
  y := (R-R0) / (R0*B);            // Zwischenvariable
  result := -x - Sqrt(Sqr(x) + y); // Temperatur berechnen
end;
 
//------------------------------------------------------------------------------
// berechnet den Widerstand bei der übergebenen Temperatur
 
 
Function CalcR(t:Double):Double;
Begin
   if t >= 0 then                    // für positive Temperaturn
    result:= R0 *(1 + A*t + B*t*t)
  else                               // negative Temperaturn
    result:= R0 *(1 + A*t + B*t*t + C*(t-100)*t*t*t);   // diese Formel möchte ich umgestellt haben nach "t"
end;
 
Grüße von Siro
Bevor ich "C" ertragen muß, nehm ich lieber Lazarus...

Benutzeravatar
af0815
Lazarusforum e. V.
Beiträge: 6198
Registriert: So 7. Jan 2007, 10:20
OS, Lazarus, FPC: FPC fixes Lazarus fixes per fpcupdeluxe (win,linux,raspi)
CPU-Target: 32Bit (64Bit)
Wohnort: Burgenland
Kontaktdaten:

Re: Temperaturumrechnung PT1000

Beitrag von af0815 »

Vielleicht im Excelsheet von hier https://www.mikrocontroller.net/topic/58407 was für dich drinnen.
Blöd kann man ruhig sein, nur zu Helfen muss man sich wissen (oder nachsehen in LazInfos/LazSnippets).

siro
Beiträge: 730
Registriert: Di 23. Aug 2016, 14:25
OS, Lazarus, FPC: Windows 11
CPU-Target: 64Bit
Wohnort: Berlin

Re: Temperaturumrechnung PT1000

Beitrag von siro »

@af0815:
Vielen Dank für den Link.
Auf jeden Fall ein sehr schönes Dokument bzw. Excel Sheet für die Berechnung ab 0 Grad.
Temperaturen unter 0 Grad werden auch berücksichtigt, wobei hier der dritte "C" Koeffizent nicht eingeht, aber ich glaube so genau brauche ich das auch nicht wirklich...
Meine Umrechnung scheint auch etwas anders zu sein, das muss ich direkt mal miteinander vergleichen.

Siro
Grüße von Siro
Bevor ich "C" ertragen muß, nehm ich lieber Lazarus...

wp_xyz
Beiträge: 4869
Registriert: Fr 8. Apr 2011, 09:01

Re: Temperaturumrechnung PT1000

Beitrag von wp_xyz »

Händisch kriegst du das nicht gelöst (obwohl ich sicher bin, dass das für die alten Meister kein Problem war). Heutzutage geht man auf wolframalpha.com und gibt in der Fragezeile den Ausdruck (ohne Anführungsstriche) "Solve t R-R0 *(1 + A*t + B*t*t + C*(t-100)*t*t*t) = 0" ein und erhält nach wenigen Sekunden die Lösung:

Code: Alles auswählen

Results:
 
t≈-0.5 sqrt(-(0.666667 B)/C + (0.264567 (2 B^3 R0^3 + 270000 C^2 R0^3 - 72 B C R0^3 + 900 B C (A R0 - R) R0^2 + 27 C (A R0 - R)^2 R0 + sqrt((2 B^3 R0^3 + 270000 C^2 R0^3 - 72 B C R0^3 + 900 B C (A R0 - R) R0^2 + 27 C (A R0 - R)^2 R0)^2 - 4 (B^2 R0^2 + 12 C R0^2 + 300 C (A R0 - R) R0)^3))^(1/3))/(C R0) + (0.419974 (R0 B^2 - 300 C R + 300 A C R0 + 12 C R0))/(C (2 B^3 R0^3 + 270000 C^2 R0^3 - 72 B C R0^3 + 900 B C (A R0 - R) R0^2 + 27 C (A R0 - R)^2 R0 + sqrt((2 B^3 R0^3 + 270000 C^2 R0^3 - 72 B C R0^3 + 900 B C (A R0 - R) R0^2 + 27 C (A R0 - R)^2 R0)^2 - 4 (B^2 R0^2 + 12 C R0^2 + 300 C (A R0 - R) R0)^3))^(1/3)) + 2500) - 0.5 sqrt(-(1.33333 B)/C - (0.264567 (2 B^3 R0^3 + 270000 C^2 R0^3 - 72 B C R0^3 + 900 B C (A R0 - R) R0^2 + 27 C (A R0 - R)^2 R0 + sqrt((2 B^3 R0^3 + 270000 C^2 R0^3 - 72 B C R0^3 + 900 B C (A R0 - R) R0^2 + 27 C (A R0 - R)^2 R0)^2 - 4 (B^2 R0^2 + 12 C R0^2 + 300 C (A R0 - R) R0)^3))^(1/3))/(C R0) - (0.419974 (R0 B^2 - 300 C R + 300 A C R0 + 12 C R0))/(C (2 B^3 R0^3 + 270000 C^2 R0^3 - 72 B C R0^3 + 900 B C (A R0 - R) R0^2 + 27 C (A R0 - R)^2 R0 + sqrt((2 B^3 R0^3 + 270000 C^2 R0^3 - 72 B C R0^3 + 900 B C (A R0 - R) R0^2 + 27 C (A R0 - R)^2 R0)^2 - 4 (B^2 R0^2 + 12 C R0^2 + 300 C (A R0 - R) R0)^3))^(1/3)) - (0.25 (-(400 B)/C - (8 (A R0 - R))/(C R0) + 1000000))/sqrt(-(0.666667 B)/C + (0.264567 (2 B^3 R0^3 + 270000 C^2 R0^3 - 72 B C R0^3 + 900 B C (A R0 - R) R0^2 + 27 C (A R0 - R)^2 R0 + sqrt((2 B^3 R0^3 + 270000 C^2 R0^3 - 72 B C R0^3 + 900 B C (A R0 - R) R0^2 + 27 C (A R0 - R)^2 R0)^2 - 4 (B^2 R0^2 + 12 C R0^2 + 300 C (A R0 - R) R0)^3))^(1/3))/(C R0) + (0.419974 (R0 B^2 - 300 C R + 300 A C R0 + 12 C R0))/(C (2 B^3 R0^3 + 270000 C^2 R0^3 - 72 B C R0^3 + 900 B C (A R0 - R) R0^2 + 27 C (A R0 - R)^2 R0 + sqrt((2 B^3 R0^3 + 270000 C^2 R0^3 - 72 B C R0^3 + 900 B C (A R0 - R) R0^2 + 27 C (A R0 - R)^2 R0)^2 - 4 (B^2 R0^2 + 12 C R0^2 + 300 C (A R0 - R) R0)^3))^(1/3)) + 2500) + 5000) + 25

So richtig programmierer-freundlich ist das auch nicht...

Ich würde stattdessen ein einfaches Nullstellenverfahren wählen und die Gleichung numerisch lösen, z.B. mit dem Newton-Verfahren (https://de.wikipedia.org/wiki/Newton-Verfahren).

Hier eine grob getestete Lösung (die Funktion "Newton" berechnet die Temperatur, die zum Widerstand R gehört:

Code: Alles auswählen

const
  R0 = 1000;
  A = 3.9083E-3;
  B = -5.775E-7;
  C = -4.183E-12;
 
// Berechnet  die Kurve R(T) für die angegeben Materialparameter R0, A, B, C
function Calc_R(T, R0, A, B, C: Double): Double;
begin
  if T >= 0 then
    Result := R0 * (1 + A*T + B*T*T)
  else
    Result := R0 * (1 + A*T + B*T*T + C*(T-100)*T*T*T);
end;
 
// Ableitung dR/dT = Steigung der R(T) Kurve
function Calc_dRdT(T, R0, A, B, C: Double): Double;
begin
  if T >= 0 then
    Result := R0 * (A + 2*B*T)
  else
    Result := R0 * (A + 2*B*T - 300*C * T*T + 4*C*T*T*T);
end;
 
function Newton(R, R0, A, B, C: Double): Double;
const
  EPS = 1E-3;   // Geforderte Genauigkeit: 0.001 Grad dürften reichen
  MAX_IT = 10// Maximalzahl der Iterationen
var
  dRdT: Double;  // Steigung der Tangente am Test-Punkt
  it: Integer;
  Rt: Double;    // Widerstand bei der Temperatur des aktuellen Iterationsschritts
  T: Double;
begin
  it := 0;
  T := 0// anfängliche Schätzung T = 0°, damit sollte es auf jeden Falle eine Lösung geben
  repeat
    Rt := Calc_R(T, R0, A, B, C);
    dRdT := Calc_dRdT(T, R0, A, B, C);
    T := T - (Rt - R)/dRdT;  // neue Lösung nach diesem Iterationsschritt
    inc(it);
  until (abs(Rt - R) < EPS) or (it > MAX_IT);
  Result := T;
end
Zuletzt geändert von wp_xyz am Di 27. Aug 2019, 23:23, insgesamt 1-mal geändert.

Benutzeravatar
six1
Beiträge: 782
Registriert: Do 1. Jul 2010, 19:01

Re: Temperaturumrechnung PT1000

Beitrag von six1 »

Respekt!
Gruß, Michael

siro
Beiträge: 730
Registriert: Di 23. Aug 2016, 14:25
OS, Lazarus, FPC: Windows 11
CPU-Target: 64Bit
Wohnort: Berlin

Re: Temperaturumrechnung PT1000

Beitrag von siro »

Ja, Respekt absolut... SUPER GEIL :wink: ich bin von den Socken, das ist doch mal eine Lösung.

wolframalpha.com kannte ich natürlich noch nicht. Sieht auch nicht wirklich soooo übersichtlich aus...

Aber ersteinmal vielen Dank wp_xyz für deine Mühe und auch den Code.

Deine Idee mit dem Newton Verfahren gefällt mir sehr gut.
Ich habe das schonmal ähnlich angewandt für eine Wurzel Berechnung, nach dem Newton/Raphson Verfahren.
Du hast das ja mein Problem sogar gleich mal in Code umgesetzt, echt supi, vielen Dank für deine Mühe und Zeit.
Ich muss mir das erstmal genauer anschauen, aber das scheint mir wirklich zur Zeit die beste Lösung sein.
Zumal "eine" Lösung und nicht "vier" wie mir zunächst suggeriert wurde... :wink:

Siro
Grüße von Siro
Bevor ich "C" ertragen muß, nehm ich lieber Lazarus...

wp_xyz
Beiträge: 4869
Registriert: Fr 8. Apr 2011, 09:01

Re: Temperaturumrechnung PT1000

Beitrag von wp_xyz »

siro hat geschrieben:Zumal "eine" Lösung und nicht "vier" wie mir zunächst suggeriert wurde... :wink:

Nein, da haben deine Freunde schon recht: Ein Polynom vierten Grades hat maximal vier Nullstellen. Aber das PT1000-Polynom ist in dem interessierenden Bereich fast linear, da gibt es nur 1 Lösung. Natürlich, wenn man bis zu Millionen Grad oder weit under den absoluten Nullpunkt geht, können die anderen Lösungen in Frage kommen, aber bis dort ist die Gleichung eh nicht mehr gültig (z.B. weil die Näherung nicht mehr stimmt, aber auch weil das Platin geschmolzen ist, oder die Tempartur nicht kleiner als -273C werden kann) - die Mathematik ist da etwas realitätsfremd...

In dem beigefügten Chart habe ich meinen Code von oben getestet: Mit dem Scrollbar den Widerstandswert (rot) einstellen, die zugehörige Temperatur wird mit dem Newton-Verfahren berechnet und an der blauenLinie angezeigt - die rote und blaue Linie müssen sich auf der berechneten R(T)-Kurve schneiden, wenn der Code für das Newton-Verfahren richtig ist.
Dateianhänge
PT1000.zip
(2.95 KiB) 279-mal heruntergeladen

siro
Beiträge: 730
Registriert: Di 23. Aug 2016, 14:25
OS, Lazarus, FPC: Windows 11
CPU-Target: 64Bit
Wohnort: Berlin

Re: Temperaturumrechnung PT1000

Beitrag von siro »

Ganz hervorragend wp_xyz,
ich habe es grad ausprobiert. Ist deckungsgleich.

Zudem ein sehr schönes Demoprojekt auch für die grafische Darstellung von Messwerten.
Hab ich gleich noch etwas dazugelernt.

Hast Dir echt nen Bier oder Pizza verdient.... Sag Becheid wenn Du mal in Berlin bist :wink:

Siro
Grüße von Siro
Bevor ich "C" ertragen muß, nehm ich lieber Lazarus...

Timm Thaler
Beiträge: 1224
Registriert: So 20. Mär 2016, 22:14
OS, Lazarus, FPC: Win7-64bit Laz1.9.0 FPC3.1.1 für Win, RPi, AVR embedded
CPU-Target: Raspberry Pi 3

Re: Temperaturumrechnung PT1000

Beitrag von Timm Thaler »

Ganz ehrlich? Da du eh nicht so genau messen kannst, nimmst du die einfachere Formel mit t und t², dann bekommst du eine quadratische Gleichung und die ist mit Tafelwerk lösbar.

Ich nehm auf dem AVR dafür eine Tabelle in 10K Schritten, zwischen denen interpoliert wird. Die Tabelle erstelle ich mit Calc.

Benutzeravatar
fliegermichl
Lazarusforum e. V.
Beiträge: 1430
Registriert: Do 9. Jun 2011, 09:42
OS, Lazarus, FPC: Lazarus Fixes FPC Stable
CPU-Target: 32/64Bit
Wohnort: Echzell

Re: Temperaturumrechnung PT1000

Beitrag von fliegermichl »

Auf der Wikipedia Seite ist u.a. zu lesen: "Das Iterationsverfahren konvergiert im günstigsten Fall asymptotisch mit quadratischer Konvergenzordnung, die Zahl der korrekten Dezimalstellen verdoppelt sich dann in jedem Schritt."

Kann mir das mal einer auf deutsch erklären?

wp_xyz
Beiträge: 4869
Registriert: Fr 8. Apr 2011, 09:01

Re: Temperaturumrechnung PT1000

Beitrag von wp_xyz »

fliegermichl hat geschrieben:Auf der Wikipedia Seite ist u.a. zu lesen: "Das Iterationsverfahren konvergiert im günstigsten Fall asymptotisch mit quadratischer Konvergenzordnung, die Zahl der korrekten Dezimalstellen verdoppelt sich dann in jedem Schritt."

Kann mir das mal einer auf deutsch erklären?

Ich denke, das heißt: Die Abweichung des gefundenen Werts wird mit jedem Iterationsschritt kleiner, und zwar stehen die Differenzen in einem quadratischen Verhältnis zueinanander. Angenommen, der berechnete Wert hat nach einem bestimmten Iterationsschritt die Differenz 0.1 zum wahren Wert. Nach dem nächsten Durchgang ist die Differenz nur noch 0.1^2 = 0.01, dann 0.01^2 = 0.0001, usw.

Hier ein WriteLn aus dem weiter oben geposteten Demoprogramm, bei dem ich die erforderliche Genauigkeit auf 1E-9 angezogen habe. Für die Berechnung der Temperatur zum Widerstand 30 Ohm (bei einem PT100, also 100 Ohm bei 0°C) gibt mir das WriteLn für jede Iteration folgende Werte aus:

Code: Alles auswählen

R = 30.000 T (direkt) = -174.601372299194050
Iteration: 0 T = 0.000000000000000
Iteration: 1 T = -179.106005168487570
Iteration: 2 T = -173.170882944500250
Iteration: 3 T = -173.157673789971260
Iteration: 4 T = -173.157673726184810

Man sieht, dass die Berechnung sehr schnell konvergiert: schon nach 5 Schritten ist die (extreme) Genauigkeitsschranke unterschritten. Als "T (direkt)" übrigens zeige ich das Ergebnis der einfachen quadratischen Formel (siro's Calc_T); diese ist für diese extreme Temperatur nur etwa 1.5° falsch. Timm Thaler's Einwand ist also nicht von der Hand zu weisen. Aber ich würde direkt die Funktion Calc_T verwenden statt einer Tabelleninterpolation. Bist du ein Fan der Logarithmentafeln?

P.S. Bei positiven Temperaturen ist die Funktion Calc_T sowieso der Newton-Methode vorzuziehen.

Antworten