1

Im working with Dymola Version 2013. I try to solve a simple mathematical problem, like:

f= x^2 -4 ;
y=1;
f=y;

f and x are defined as Real. The solution is 2.36. but I need to calculate with both solutions. So 2.36 and -2.36! In my problem f is a polynomial like ax^3 + bx^2 +cx +d and y is linear. y = ax +b

How can I get all solutions of this problem? x has no explicit value. x often has at least two solutions. Should x be a vector? in this case I got problems with the dimension of the equation... Can someone help me?

  • Did you already take a look at `Modelica.Math.Vectors.Utilities.roots`? It is a function to calculate the roots of a polynomial. If it helps, you should write a short example as a self-answer here. – matth Apr 14 '13 at 15:40

2 Answers2

1

If I understand you correctly, your goal is to use Modelica to find all the roots of a (higher than second order) polynomial. I'm afraid that just isn't what Modelica is intended for. For a given non-linear equations, the simulation of a Modelica model will use (at most) one root of the non-linear equation. If you want to find all the roots, you'll have to find a way to factor the polynomial yourself. In your case, you are only dealing with a cubic polynomial so you should research algorithms for factoring cubic polynomials. You could then write such an algorithm as a Modelica function.

Michael Tiller
  • 9,291
  • 3
  • 26
  • 41
1

As I understood your question you have two polynomials and want to find all points where they are equal.
Here is a function that does that using Modelica.Math.Vectors.Utilities.roots:
First, you give the two polynomials poly1 and poly2. Finding poly1=poly2is identical to finding poly1-poly2=0, so I define a third polynomial polyDiff = polyLong-polyShort and then hand over that polynomial to Modelica.Math.Vectors.Utilities.roots. It will return all roots, even complex ones.

function polyIntersect
  input Real[:] poly1={3,2,1,0};
  input Real[:] poly2={8,7};
  output Real[:,2] intersect;

protected 
  Integer nPoly1 = size(poly1,1);
  Integer nPoly2 = size(poly2,1);
  Integer nPolyShort = min(nPoly1, nPoly2);
  Integer nPolyLong = max(nPoly1, nPoly2);
  Real[nPolyShort] polyShort;
  Real[nPolyLong] polyLong;
  Real[nPolyLong] polyDiff;

algorithm 
  if (nPoly1<nPoly2) then 
    polyShort := poly1;
    polyLong := poly2;
  else 
    polyShort := poly2;
    polyLong := poly1;
  end if;

  polyDiff := polyLong;
  for i in 0:nPolyShort-1 loop
    polyDiff[nPolyLong-i] := polyLong[nPolyLong-i] - polyShort[nPolyShort-i];
  end for;

  intersect := Modelica.Math.Vectors.Utilities.roots(polyDiff);

end polyIntersect;

The above code is also available here: https://gist.github.com/thorade/5388205

matth
  • 2,568
  • 2
  • 21
  • 41