1. Limited time only! Sign up for a free 30min personal tutor trial with Chegg Tutors
    Dismiss Notice
Dismiss Notice
Join Physics Forums Today!
The friendliest, high quality science and math community on the planet! Everyone who loves science is here!

C++: Method of Minimum Squares used for System Solving

  1. Apr 9, 2017 #1
    1. The problem statement, all variables and given/known data

    Okay, this one is a bit big, and I'm attempting a translation, so I'm gonna post a TL;DR version at the bottom just to be sure. Anyway, here goes:

    While doing an experiment, we write down the values of a physical/natural size y, for various values of a physical/natural size x. Let's say n is the amount of (x,y) pairs which are the result of how many times we did this experiment. We know that these two sizes are connected through the following formula: y = f(x) = a*x + b
    Sadly, it's easy to have lots of faults with such a way of writting our measurments. So, we need to find a linear relation between x & y (so essentially the values of a & b), which is as close as it can to (xi,yi).
    So we need to find the formula f(x) = a*x + b that minimizes the average square fault E = √(1/n ** ∑i=1n(f(xi) - yi)2)
    This is equilavent to this system:

    n*b + (∑i=1nxi)*a = ∑i=1nyi

    (∑i=1nxi)*b + (∑i=1nxi2)*a = ∑i=1nxi*yi

    You have to write a program in C++ & Fortran that get close to the linear relation of x & y, by using the pairs of (xi,yi). The amount of measurements (n) and the pairs themselves will be entered through the keyboard, by the user.

    Use the Program from Exercise 2 to solve the system.


    TLD;DR: Create a program that takes n amount of (x,y) pairs, entered through the keyboard, and solves the above equation system, with a & b being the unknown quantities.

    PS: That floating n over at the ∑ is just the "limit", it's not Σ to the nth power.

    2. Relevant equations


    3. The attempt at a solution

    Okay, as far as this program goes, I have this:

    Code (C):


    #include <iostream>
    #include <cmath>
    #include <cstdlib>
    using namespace std;

    int main()
    {
     int i, j;
     int n;
     double coor[50][2];

     cout << "How many measurments did you take ? ";
     cin >> n;

     for(i=1; i<=n; i++)
      for(j=1; j<=n; j++)
       {
        cout << "Enter your coordinates for coor(" << i << ") (" << j << "): ";
        cin >> coor[i][j];
       }
       
       return 0;
    }

     
    Yeah, it's nothing thus far, but I can't "solve" the input problem yet, so I'm stumped.

    The program that solves the system, from a previous excercise:

    Code (C):


    #include <iostream>
    using namespace std;

    int main()
    {
     double a11, a12, a21, a22;
     double b1, b2;
     double x, y, D, Dx, Dy;

     cout << "Give the values of a11, a12, a21 & a22: " << endl;
     cin >> a11 >> a12 >> a21 >> a22;
     cout << endl;

     cout << "Give the values of b1 & b2: " << endl;
     cin >> b1 >> b2;
     cout << endl;

     D = (a11*a22 - a12*a21);
     Dx = (b1*a22 - b2*a12);
     Dy = (b2*a11 - b1*a21);

    if (D != 0)
    {
     x = Dx/D;
     y = Dy/D;
     cout << "For those values, x & y are: " << x << " & " << y << endl;
    }

    if (D = 0)
    {
     if (Dx = 0)
     {
      if (Dy = 0)
      {
       cout << "Then the equations/system have infinite solutions." << endl;
      }
     }
    }

    if (D = 0)
    {
     if (Dx != 0)
     {
      if (Dy != 0)
       {
        cout << "The equations/system are unsolvable." << endl;
       }
      }
    }

    return 0;
    }
       
     
    It's a tad different, but the basic idea is the same.

    Anyway, I'd greatly appreciate any kind of help!
     
  2. jcsd
  3. Apr 9, 2017 #2

    BvU

    User Avatar
    Science Advisor
    Homework Helper
    Gold Member

    Ah, you mean $$\sum_{i=1}^n$$ In PF that's written as $$\sum_{i=1}^n$$ !

    If you enter e.g. n = 10 you want to go through the
    Code (C):
      {
        cout << "Enter your coordinates for coor(" << i << ") (" << j << "): ";
        cin >> coor[i][j];
       }
     
    loop how many times ? As presented, i runs from 1 to n and for each i, j runs from 1 to n too, so that is n2 times !
     
  4. Apr 9, 2017 #3
    Thanks! I didn't know about that, so I just improvised. I'll keep it in mind for next time.

    Yeah, I noticed the problem when I run it, but I'm out of ideas. What should happen, would be this:

    E.G. n = 3:
    Enter your coordinates for coor[1][1] = ...
    Enter your coordinates for coor[1][2] = ...
    Enter your coordinates for coor[2][1] = ...
    Enter your coordinates for coor[2][2] = ...
    Enter your coordinates for coor[3][1] = ...
    Enter your coordinates for coor[3][2] = ...

    But, the result is obviously not that. I figure the problem is somewhere in the j <= ... part. If I put in i, it disregards one half of the array. if I put in n, well, it creates many more "columns". If I put in 2 (since I want my array to have one spot for x, and one for y, essentially), it stops before I'm able to put in all of my "coordinates"/data for each measurement.
     
  5. Apr 9, 2017 #4
    arrays start at 0, so you should have j=0; j<=1 in the for loop.
     
  6. Apr 9, 2017 #5
    I know, but I was under the impression that if I declared it, I could start at 1 or whatever other number. I've not really worked with arrays before, but I figured since I created room for a big array, even if the input for n was something like, say, 25, there'd be plenty of room to start from 1 and go coor[1][1] -> coo[1][2] etc. Correct me if I'm wrong though, like I said, I've had minimal hands-on with arrays.
     
  7. Apr 9, 2017 #6

    BvU

    User Avatar
    Science Advisor
    Homework Helper
    Gold Member

    So you only want to go through n times. Easiest solution: one array for x and one array for y !
     
  8. Apr 9, 2017 #7
    Yeah, that's actually a nifty way around it. I figured I'd get some work done with multidimensional arrays to get the hang of them, but I suppose I should familiarize myself with the simple ones first.

    Okay, it's midnight here, so I'm gonna try the program tomorrow and get back to the thread. Thanks for the help so far everyone!
     
  9. Apr 10, 2017 #8
    I solved it on C++, but I have a problem with the Fortran Version, about halfway through.

    First up, the working C++ version:

    Code (C):


    #include <iostream>
    #include <cmath>
    #include <cstdlib>
    using namespace std;

    int main()
    {
     int i, j;
     int n;
     double coorx[50];
     double coory[50];
     double sumx, sumy, sumxs, sumxy;
     double a, b;
     double D, Da, Db;

     cout << "How many measurments did you take ? ";
     cin >> n;
     cout << endl;

     sumx = 0;
     sumy = 0;
     sumxs = 0;
     sumxy = 0;

     for(i=1; i<=n; i++)
      {
       cout << "Enter the x coordinate for Measurement No " << i << " in coorx(" << i << "): ";
       cin >> coorx[i];
       sumx += coorx[i];
       sumxs += (coorx[i]*coorx[i]);
      }
     
      cout << endl << "The sum of x coordinates is: " << sumx << endl << "And the sum of squared coordinates x is: " << sumxs << endl << endl;
     
      for(j=1; j<=n; j++)
      {
       cout << "Enter the y coordinate for Measurement No " << j << " in coory(" << j << "): ";
       cin >> coory[j];
       sumy += coory[j];
      }
     
      cout << endl << "The sum of y coordinates is: " << sumy << endl << endl;
     
      for(i=1, j=1; i<=n, j<=n; i++, j++)
      {
       sumxy += coorx[i]*coory[j];
      }
       
      cout << "The sum of the products between coordinates x and y for each measurement is: " << sumxy << endl << endl;
     
     D = (n*sumxs - sumx*sumx);
     Db = (sumy*sumxs - sumxy*sumx);
     Da = (sumxy*n - sumy*sumx);

    if (D != 0)
    {
     b = Db/D;
     a = Da/D;
     cout << "For those values, a & b are: " << a << " & " << b << endl << endl;
    }

    if (D = 0)
    {
     if (Da = 0)
     {
      if (Db = 0)
      {
       cout << "Then the equations/system have infinite solutions." << endl << endl;
      }
     }
    }

    if (D = 0)
    {
     if (Da != 0)
     {
      if (Db != 0)
       {
        cout << "The equations/system are unsolvable." << endl << endl;
       }
      }
    }

    return 0;
    }
       
     
    And now here's my Fortran one:

    Code (Fortran):


    PROGRAM SquarMeth

    REAL, DIMENSION(50) :: coorx
    REAL, DIMENSION(50) :: coory

    INTEGER :: i, j
    INTEGER :: n

    REAL :: sumx, sumy, sumxs, sumxy
    REAL :: D, Da, Db
    REAL :: a, b

    PRINT *, 'How many measurements did you take ?'
    READ *, n
    PRINT *

    sumx = 0
    sumy = 0
    sumxs = 0
    sumxy = 0

    DO i=1,n
    PRINT *,'Enter the x coordinate for Measurement No', i ,'in coorx(',i,'):'
    READ *, coorx(i)
    sumx = sumx + coorx(i)
    sumxs = sumxs + (coorx(i)*coorx(i))
    ENDDO

    PRINT *
    PRINT *,'The sum of x coordinates is:',sumx
    PRINT *
    PRINT *,'And the sum of squared coordinates x is:',sumxs
    PRINT *
    PRINT *

    DO j=1,n
    PRINT *,'Enter the y coordinate for Measurement No',j,'in coory(',j,'):'
    READ *, coory(j)
    sumy = sumy + coory(j)
    ENDDO

    PRINT *
    PRINT *,'The sum of y coordinates is:',sumy
    PRINT *
    PRINT *

    DO i=1,n
    DO j=1,n
    sumxy = sumxy + coorx(i)*coory(j)
    ENDDO
    ENDDO

    PRINT *
    PRINT *,'The sum of the products between coordinates x and y for each measurement is:',sumxy
    PRINT *
    PRINT *

    D = (n*sumxs - sumx*sumx)
    Db = (sumy*sumxs - sumxy*sumx)
    Da = (sumxy*n - sumy*sumx)

    IF (D /= 0) THEN
    b = Db / D
    a = Da / D
    PRINT *, 'For these values, b & a are:', b,'&',a
    END IF

    IF (D == 0) THEN
     IF (Db == 0) THEN
      IF (Da == 0) THEN
       PRINT *, 'The equations/system have infinite solutions.'
      END IF
     END IF
    END IF

    IF (D == 0) THEN
     IF (Db /= 0) THEN
      IF (Da /= 0) THEN
       PRINT *, 'The equations/system are unsolvable.'
      END IF
     END IF
    END IF

    END PROGRAM SquarMeth

     
    The problem arises when it comes to the program having to compute sumxy, which skews the results of the rest of the "variables" and messes up the program.

    Any help is appreciated!
     
  10. Apr 10, 2017 #9

    BvU

    User Avatar
    Science Advisor
    Homework Helper
    Gold Member

    Same as before: you don't want ##\displaystyle{\sum_i\sum_j\operatorname {coorx}(i) \operatorname {coory}(j) \ \ }## (n2 terms), but you want ##\ \ \displaystyle\sum_i\operatorname {coorx}(i) \operatorname {coory}(i)\ \ ## (only n terms)

    PS wouldn't it be easier for the user to enter x,y pairs instead of first a whole x array and then a whole y array ? i.e. use
    Code (Fortran):
    DO i=1,n
    PRINT *,'Enter the x and y coordinates for Measurement No', i ,'):'
    READ *, coorx(i), coory(i)

    sumx = sumx + coorx(i)
    sumxs = sumxs + (coorx(i)*coorx(i))
    sumy = ...
    sumys ...
    sumxy ...
    ENDDO
    [edit]This way, the C version would become easier to read too
     
    Last edited: Apr 10, 2017
  11. Apr 10, 2017 #10
    Ah yeah, I got it (I tried it and it worked). Thanks!

    PS: Yeah, it's a bit "drawn out", but I like things being more "detailed", with more cout(s) and whatnot. I know it's not recommended but considering these are home-programs that nobody really checks, I figure it's okay to do things a bit "my way" until I get the hang of it.

    Again, thanks a ton for the help!
     
  12. Apr 10, 2017 #11

    Mark44

    Staff: Mentor

    Your C++ code has some glaring errors that no one mentioned. Neither of the two outer if blocks above will ever execute. The expression D = 0 is not comparing the value of D with zero -- it is assigning 0 to D. To compare two things for equality, you need to use the equality operator, ==, not the assignment operator, =.

    In your code, if D happens to be, say 3 before entering the first if block, if (D = 0) sets D to 0. Since 0 is considered false, control flows to the second outer if block. Since D has already been set to 0, that test expression is still considered false, so control flows to whatever statements are after the second outer if block, which in your code would be the return 0; statement.
     
  13. Apr 10, 2017 #12
    Oh darn, you're right. It compiled and worked with certain values, that I didn't catch that part. Thanks a ton!
     
  14. Apr 10, 2017 #13

    Mark44

    Staff: Mentor

    The compiler can catch only syntax errors -- errors for which you have incorrect language usage. It can't catch errors in your logic (semantic errors) such as the ones I pointed out.

    Also, your whole if statement contains some duplication that could be eliminated.

    Code (C):
    if (D != 0)
    {
        b = Db/D;
        a = Da/D;
        cout << "For those values, a & b are: " << a << " & " << b << endl << endl;
    }
    else     // D == 0
    {
        if (Da == 0 && Db == 0)
        {
           cout << "Then the equations/system have infinite solutions." << endl << endl;
        }
        if (Da != 0 && Db != 0)
        {
           cout << "The equations/system are unsolvable." << endl << endl;
        }
    }
    I have merely cleaned up your logic. I haven't looked closely enough at the whole program to say that your logic takes care of all possibilities. In particular, the above control structure doesn't deal with the cases where Da == 0 and Db != 0, or when Da != 0 and Db == 0. In the else clause above, either both Da and Db are zero, or they are both nonzero.
     
  15. Apr 11, 2017 #14
    The "equation system" part was a remnant from a previous program that we were to use, so I did a copy/paste and didn't really look into it. Thanks for the cleaning up and the piece of advice!
     
Know someone interested in this topic? Share this thread via Reddit, Google+, Twitter, or Facebook

Have something to add?
Draft saved Draft deleted



Similar Discussions: C++: Method of Minimum Squares used for System Solving
  1. Newton's method in c++ (Replies: 8)

Loading...