pow() with a negative base
MPD78
10-19-2010, 10:36 AM
Hello all,
I have an application that requires using a negative base with the pow() function. This application does not allow for simply passing the abs() of the base and then multiplying the returned pow() result by -1.
What is a work around for this?
Thanks
Matt
davekw7x
10-19-2010, 05:21 PM
We are talking about real numbers, not complex numbers, right?
Well...
The "schoolboy" definition of a number raised to a power goes something like this:
n^m = \textrm{ ``n multiplied by itself m times''}
That works OK if n and m are both positive integers. It even works if n is a positive number (not necessarily an integer) and m is a positive integer. It even works if n is a negative integer and m is a non-negative integer. It also "turns out" that if m is a negative integer, you can just take 1 divided by (n multiplied by itself abs(m) times) to be a reasonable result.
But how would you apply that definition to something like (e raised to the pi power)? How the heck can you multiply 2.718... by itself 3.141... times?
So, welcome to the "real" world where a number raised to a power is defined by something like the following for positive values of a
a^b = \exp(b \log a)
The variables a and b don't necessarily have to be integers.
However...
If a is a negative (or zero), then that expression is not defined, since the logarithm of a number is defined only for positive numbers. The variable b can be any kind of real number. (Integer, non-integer, positive, negative, zero.)
Note that if m is a negative integer and n is a positive integer, you can use the first definition to get an answer even though the "real" definition doesn't hold.
Anyhow...
The standard math function double pow(double, double) takes some of these special cases into account:
#include <iostream>
#include <iomanip>
#include <cmath>
using namespace std;
int main()
{
cout << fixed << setprecision(3);
for (int i = -5; i < 6; i++) {
for (int x = 2; x <= 3; x++) {
cout << setw(6) << i << " to the power " << x
<< " = " << setw(8) << pow((double)i, (double)x) << endl;
}
}
cout << endl;
for (double i = -5.1; i < 6; i++) {
double x = 3.14;
cout << setw(6) << i << " to the power " << x
<< " = " << setw(8) << pow(i, x) << endl;
}
cout << endl;
for (int i = -5; i < 6; i++) {
double x = 3.14;
cout << setw(6) << i << " to the power " << x
<< " = " << setw(8) << pow((double)i, x) << endl;
}
cout << endl;
return 0;
}
Output (GNU g++ compiler). Note that the "nan" values in the output listing stand for "not a number," meaning that the pow() function knows that you can't do that, and, in fact, doesn't even try. (Older compilers might not have been so careful, and might have given very large or very small numbers with no indication that they were erroneous, but the treatment of "nan" things is part of current C and C++ standards.)
-5 to the power 2 = 25.000
-5 to the power 3 = -125.000
-4 to the power 2 = 16.000
-4 to the power 3 = -64.000
-3 to the power 2 = 9.000
-3 to the power 3 = -27.000
-2 to the power 2 = 4.000
-2 to the power 3 = -8.000
-1 to the power 2 = 1.000
-1 to the power 3 = -1.000
0 to the power 2 = 0.000
0 to the power 3 = 0.000
1 to the power 2 = 1.000
1 to the power 3 = 1.000
2 to the power 2 = 4.000
2 to the power 3 = 8.000
3 to the power 2 = 9.000
3 to the power 3 = 27.000
4 to the power 2 = 16.000
4 to the power 3 = 64.000
5 to the power 2 = 25.000
5 to the power 3 = 125.000
-5.100 to the power 3.140 = nan
-4.100 to the power 3.140 = nan
-3.100 to the power 3.140 = nan
-2.100 to the power 3.140 = nan
-1.100 to the power 3.140 = nan
-0.100 to the power 3.140 = nan
0.900 to the power 3.140 = 0.718
1.900 to the power 3.140 = 7.504
2.900 to the power 3.140 = 28.309
3.900 to the power 3.140 = 71.770
4.900 to the power 3.140 = 146.966
5.900 to the power 3.140 = 263.315
-5 to the power 3.140 = nan
-4 to the power 3.140 = nan
-3 to the power 3.140 = nan
-2 to the power 3.140 = nan
-1 to the power 3.140 = nan
0 to the power 3.140 = 0.000
1 to the power 3.140 = 1.000
2 to the power 3.140 = 8.815
3 to the power 3.140 = 31.489
4 to the power 3.140 = 77.708
5 to the power 3.140 = 156.591
(Left as an exercise for the student: Try the same kind of cases with negative exponent values.)
Bottom line: What is the nature of the numbers in your expression? Integers? Non-integers? What?
Regards,
Dave
Footnote:
My first set of numbers illustrates something that we observed with integers raised to integer powers: A negative number raised to an odd power is a negative integer and a negative number raised to an even power is an positive integer. The pow() function takes these special cases into account rather than choking by trying to take the logarithm of a negative number. Clever, yes???
MPD78
10-20-2010, 07:17 AM
We are talking about real numbers, not complex numbers, right?
Yes, we are talking about real numbers.
Here is the output from the code you provided above.
-5 to the power 2 = 25.000
-5 to the power 3 = -125.000
-4 to the power 2 = 16.000
-4 to the power 3 = -64.000
-3 to the power 2 = 9.000
-3 to the power 3 = -27.000
-2 to the power 2 = 4.000
-2 to the power 3 = -8.000
-1 to the power 2 = 1.000
-1 to the power 3 = -1.000
0 to the power 2 = 0.000
0 to the power 3 = 0.000
1 to the power 2 = 1.000
1 to the power 3 = 1.000
2 to the power 2 = 4.000
2 to the power 3 = 8.000
3 to the power 2 = 9.000
3 to the power 3 = 27.000
4 to the power 2 = 16.000
4 to the power 3 = 64.000
5 to the power 2 = 25.000
5 to the power 3 = 125.000
-5.100 to the power 3.140 = -1.#IO
-4.100 to the power 3.140 = -1.#IO
-3.100 to the power 3.140 = -1.#IO
-2.100 to the power 3.140 = -1.#IO
-1.100 to the power 3.140 = -1.#IO
-0.100 to the power 3.140 = -1.#IO
0.900 to the power 3.140 = 0.718
1.900 to the power 3.140 = 7.504
2.900 to the power 3.140 = 28.309
3.900 to the power 3.140 = 71.770
4.900 to the power 3.140 = 146.966
5.900 to the power 3.140 = 263.315
-5 to the power 3.140 = -1.#IO
-4 to the power 3.140 = -1.#IO
-3 to the power 3.140 = -1.#IO
-2 to the power 3.140 = -1.#IO
-1 to the power 3.140 = -1.#IO
0 to the power 3.140 = 0.000
1 to the power 3.140 = 1.000
2 to the power 3.140 = 8.815
3 to the power 3.140 = 31.489
4 to the power 3.140 = 77.708
5 to the power 3.140 = 156.591
Here is the code that I have minus the details. (If you want to compile this, let me know and I can give the details.)
// Define and load the matrix
Matrix<double,2> CHMNH2O(3,3);
CHMNH2O[0][0] = -2.2118;
CHMNH2O[0][1] = -1.1987;
CHMNH2O[0][2] = 0.035596;
CHMNH2O[1][0] = 0.85667;
CHMNH2O[1][1] = 0.93048;
CHMNH2O[1][2] = -0.14391;
CHMNH2O[2][0] = -0.10838;
CHMNH2O[2][1] = -0.17156;
CHMNH2O[2][2] = 0.045915;
for(int i=0; i<3; i++){
for(int j = 0; j<3; j++){
jj = CHMNH2O[i][j];
EPSO1 = pow(t,jj);
EPSO2 = pow(PALL,jj); // PALL is negative.
EPSO += jj*EPSO1*EPSO2;
}
}
As you can see, I have a case where NaN is being returned.
Thanks
Matt
davekw7x
10-20-2010, 09:22 AM
...
negative base with the pow() function.
What is a work around for this?
...we are talking about real numbers.....
.
.
-5.100 to the power 3.140 = -1.#IO
...
jj = CHMNH2O[i][j];//<--- not an int, right?
...
EPSO2 = pow(PALL,jj); // PALL is negative.
...
As you can see, I have a case where NaN is being returned.
Hmmm...
If you have an application that requires a negative number to be raised to a non-integer power and give a real number, well, you are SOL (Somewhat Out of Luck). Needs a little rethink.
On the other hand if we were talking about complex numbers:
#include <iostream>
#include <iomanip>
#include <cmath>
#include <complex>
using namespace std;
int main()
{
cout << scientific << showpos;
for (double i = -5; i < 6; i++) {
complex<double> ci(i);
for (double x = 2; x <= 3; x++) {
cout << ci << " to the power " << x
<< " = " << pow(ci, x) << endl;
}
}
cout << endl;
for (double i = -5.1; i < 6; i++) {
complex<double> ci(i);
double x = 3.14;
cout << ci << " to the power " << x
<< " = " << pow(ci, x) << endl;
}
cout << endl;
for (double i = -5; i < 6; i++) {
complex<double> ci(i);
double x = 3.14;
cout << ci << " to the power " << x
<< " = " << pow(ci, x) << endl;
}
cout << endl;
return 0;
}
Output:
(-5.000000e+00,+0.000000e+00) to the power +2.000000e+00 = (+2.500000e+01,-6.123032e-15)
(-5.000000e+00,+0.000000e+00) to the power +3.000000e+00 = (-1.250000e+02,+4.592274e-14)
(-4.000000e+00,+0.000000e+00) to the power +2.000000e+00 = (+1.600000e+01,-3.918740e-15)
(-4.000000e+00,+0.000000e+00) to the power +3.000000e+00 = (-6.400000e+01,+2.351244e-14)
(-3.000000e+00,+0.000000e+00) to the power +2.000000e+00 = (+9.000000e+00,-2.204291e-15)
(-3.000000e+00,+0.000000e+00) to the power +3.000000e+00 = (-2.700000e+01,+9.919311e-15)
(-2.000000e+00,+0.000000e+00) to the power +2.000000e+00 = (+4.000000e+00,-9.796851e-16)
(-2.000000e+00,+0.000000e+00) to the power +3.000000e+00 = (-8.000000e+00,+2.939055e-15)
(-1.000000e+00,+0.000000e+00) to the power +2.000000e+00 = (+1.000000e+00,-2.449213e-16)
(-1.000000e+00,+0.000000e+00) to the power +3.000000e+00 = (-1.000000e+00,+3.673819e-16)
(+0.000000e+00,+0.000000e+00) to the power +2.000000e+00 = (+0.000000e+00,+0.000000e+00)
(+0.000000e+00,+0.000000e+00) to the power +3.000000e+00 = (+0.000000e+00,+0.000000e+00)
(+1.000000e+00,+0.000000e+00) to the power +2.000000e+00 = (+1.000000e+00,+0.000000e+00)
(+1.000000e+00,+0.000000e+00) to the power +3.000000e+00 = (+1.000000e+00,+0.000000e+00)
(+2.000000e+00,+0.000000e+00) to the power +2.000000e+00 = (+4.000000e+00,+0.000000e+00)
(+2.000000e+00,+0.000000e+00) to the power +3.000000e+00 = (+8.000000e+00,+0.000000e+00)
(+3.000000e+00,+0.000000e+00) to the power +2.000000e+00 = (+9.000000e+00,+0.000000e+00)
(+3.000000e+00,+0.000000e+00) to the power +3.000000e+00 = (+2.700000e+01,+0.000000e+00)
(+4.000000e+00,+0.000000e+00) to the power +2.000000e+00 = (+1.600000e+01,+0.000000e+00)
(+4.000000e+00,+0.000000e+00) to the power +3.000000e+00 = (+6.400000e+01,+0.000000e+00)
(+5.000000e+00,+0.000000e+00) to the power +2.000000e+00 = (+2.500000e+01,+0.000000e+00)
(+5.000000e+00,+0.000000e+00) to the power +3.000000e+00 = (+1.250000e+02,+0.000000e+00)
(-5.100000e+00,+0.000000e+00) to the power +3.140000e+00 = (-1.507773e+02,-7.095041e+01)
(-4.100000e+00,+0.000000e+00) to the power +3.140000e+00 = (-7.598133e+01,-3.575410e+01)
(-3.100000e+00,+0.000000e+00) to the power +3.140000e+00 = (-3.158212e+01,-1.486142e+01)
(-2.100000e+00,+0.000000e+00) to the power +3.140000e+00 = (-9.296814e+00,-4.374749e+00)
(-1.100000e+00,+0.000000e+00) to the power +3.140000e+00 = (-1.220502e+00,-5.743248e-01)
(-1.000000e-01,+0.000000e+00) to the power +3.140000e+00 = (-6.554893e-04,-3.084498e-04)
(+9.000000e-01,+0.000000e+00) to the power +3.140000e+00 = (+7.183258e-01,+0.000000e+00)
(+1.900000e+00,+0.000000e+00) to the power +3.140000e+00 = (+7.503887e+00,+0.000000e+00)
(+2.900000e+00,+0.000000e+00) to the power +3.140000e+00 = (+2.830934e+01,+0.000000e+00)
(+3.900000e+00,+0.000000e+00) to the power +3.140000e+00 = (+7.176999e+01,+0.000000e+00)
(+4.900000e+00,+0.000000e+00) to the power +3.140000e+00 = (+1.469656e+02,+0.000000e+00)
(+5.900000e+00,+0.000000e+00) to the power +3.140000e+00 = (+2.633148e+02,+0.000000e+00)
(-5.000000e+00,+0.000000e+00) to the power +3.140000e+00 = (-1.416875e+02,-6.667305e+01)
(-4.000000e+00,+0.000000e+00) to the power +3.140000e+00 = (-7.031273e+01,-3.308666e+01)
(-3.000000e+00,+0.000000e+00) to the power +3.140000e+00 = (-2.849222e+01,-1.340742e+01)
(-2.000000e+00,+0.000000e+00) to the power +3.140000e+00 = (-7.976268e+00,-3.753347e+00)
(-1.000000e+00,+0.000000e+00) to the power +3.140000e+00 = (-9.048271e-01,-4.257793e-01)
(+0.000000e+00,+0.000000e+00) to the power +3.140000e+00 = (+0.000000e+00,+0.000000e+00)
(+1.000000e+00,+0.000000e+00) to the power +3.140000e+00 = (+1.000000e+00,+0.000000e+00)
(+2.000000e+00,+0.000000e+00) to the power +3.140000e+00 = (+8.815241e+00,+0.000000e+00)
(+3.000000e+00,+0.000000e+00) to the power +3.140000e+00 = (+3.148914e+01,+0.000000e+00)
(+4.000000e+00,+0.000000e+00) to the power +3.140000e+00 = (+7.770847e+01,+0.000000e+00)
(+5.000000e+00,+0.000000e+00) to the power +3.140000e+00 = (+1.565906e+02,+0.000000e+00)
Note that raising complex numbers to floating point powers might result in roundoff, so that (-5.0, 0.0) to the power 2.0 gives a complex number whose imaginary part is not exactly equal to zero. That's the way it goes with floating point arithmetic.
Bottom line: what would you like -5.0 to the power 3.14 to be? What kind of "workaround" would you like?
Regards,
Dave
MPD78
10-21-2010, 07:15 AM
The work around I am looking for is to be able to solve the equation that has this negative base.
Any ideas?
Thanks for you help.
Matt
davekw7x
10-21-2010, 09:36 AM
...solve the equation that has this negative base....
What equation?
If the solution does stuff like raising a negative number to a non-integer power, it will involve complex numbers, as I tried to show.
Regards,
Dave
MPD78
10-21-2010, 10:34 AM
Sorry, Dave, I forgot to attach the pdf with the equation.
Here is my code.
Note: Matix.h is not included. (I can send it to you if you like.)
#include "C:\NR Headers\nr3.h"
#include "matrix.h".
using namespace Numeric_lib;
struct GAS_RAD{
// Object used to create all the needed gas radiation calculations
vector<double> vol;
vector<int> num;
double temp;
double press_tot;
double xl;
GAS_RAD(const vector<double> vvol, const vector<int> nnum, double ttemp, double ppress_tot, double xxl)
: vol(vvol), num(nnum), temp(ttemp), press_tot(ppress_tot), xl(xxl) {} // Constructor
double epsh2o();
double epsco2();
double epstot();
};
double GAS_RAD::epsh2o(){
// Function to calculate total emissivity of water vapor (steam)
double volume = 0;
// Search the gas composition for H2O and its volume
for(int i = 0; i<vol.size(); i++){
for(int j=0; j<vol.size(); j++){
if(num[i]==3) volume = vol[j];
}
}
// Convert temperature from F to K
temp = temp*0.55555 + 290.93;
// Convert pressure from psig to Bar
double press_tot_bar = press_tot/14.50;
// Calculate the partial pressure
double part_press_h2o = volume*press_tot_bar;
// Define and load the matrix
Matrix<double,2> CHMNH2O(3,3);
CHMNH2O[0][0] = -2.2118;
CHMNH2O[0][1] = -1.1987;
CHMNH2O[0][2] = 0.035596;
CHMNH2O[1][0] = 0.85667;
CHMNH2O[1][1] = 0.93048;
CHMNH2O[1][2] = -0.14391;
CHMNH2O[2][0] = -0.10838;
CHMNH2O[2][1] = -0.17156;
CHMNH2O[2][2] = 0.045915;
double t = temp/1000;
if(t < 0.75) double a = 2.144;
double a = 1.888 - 2.053*log10(t);
double b = 1.10/(pow(t,1.4));
double c = 0.5;
double PALM = 13.2*(pow(t,2));
double PE = press_tot+2.56*part_press_h2o/sqrt(t);
double PAL = part_press_h2o*xl;
double PALL = log10(PAL);
double EPSH2O = 0; // initialize to zero
double EPSO1 = 0;
double EPSO2 = 0;
double EPSO = 0;
double jj = 0;
for(int i=0; i<3; i++){
for(int j = 0; j<3; j++){
jj = CHMNH2O[i][j];
EPSO1 = pow(t,jj);
EPSO2 = pow(PALL,jj); // PALL is negative.
EPSO += jj*EPSO1*EPSO2;
}
}
EPSO = exp(EPSO);
double EPSR = 1.0 - (a-1.0)*(1.0-PE)/(a+b-1.0+PE)*exp(-c*pow((log10(PALM)-PALL),2));
return EPSO*EPSR;
}
See attachments.
Thanks
Matt.
davekw7x
10-21-2010, 01:28 PM
...pdf with the equation.
Now it's there, but not real clear.
However...
As far as I can tell: Under the double summation, there is something in parentheses raised to the power j, and there is something else in parentheses raised to the power i, where i and j are the (integer) summation index variables. Or am I missing something?
The standard C++ pow() function can handle negative floating point numbers raised to an integer power. (See Footnote.)
Here is my code.
I didn't look real close to try to figure it out, but since you are raising a negative number to a non-integer power, I don't believe it is the proper way to implement the evaluation of the equation in the paper (or at least not the thing that I think I see in your attachment).
Regards,
Dave
Footnote:
I guess that I didn't include an example of negative floating point numbers raised to an integer power in my previous posts, so here goes:
#include <iostream>
#include <iomanip>
#include <cmath>
using namespace std;
int main()
{
cout << fixed << setprecision(3);
for (double a = -3.1; a < 3; a++) {
for (int b = 0; b < 3; b++) {
cout << setw(6) << a << " to the power " << b
<< " = " << setw(8) << pow(a, b) << endl;
}
cout << endl;
}
cout << endl;
return 0;
}
Output (GNU g++, Microsoft cl.exe, Borland bcc32.exe):
-3.100 to the power 0 = 1.000
-3.100 to the power 1 = -3.100
-3.100 to the power 2 = 9.610
-2.100 to the power 0 = 1.000
-2.100 to the power 1 = -2.100
-2.100 to the power 2 = 4.410
-1.100 to the power 0 = 1.000
-1.100 to the power 1 = -1.100
-1.100 to the power 2 = 1.210
-0.100 to the power 0 = 1.000
-0.100 to the power 1 = -0.100
-0.100 to the power 2 = 0.010
0.900 to the power 0 = 1.000
0.900 to the power 1 = 0.900
0.900 to the power 2 = 0.810
1.900 to the power 0 = 1.000
1.900 to the power 1 = 1.900
1.900 to the power 2 = 3.610
2.900 to the power 0 = 1.000
2.900 to the power 1 = 2.900
2.900 to the power 2 = 8.410
MPD78
10-21-2010, 01:52 PM
As far as I can tell: Under the double summation, there is something in parentheses raised to the power j, and there is something else in parentheses raised to the power i, where i and j are the (integer) summation index variables. Or am I missing something?
That is correct. Sorry for the poor resolution.
Thanks
Matt
davekw7x
10-21-2010, 05:08 PM
That is correct....
So, the point that I was trying to make is this:
Since the mathematical expression involves raising (sometimes-negative) stuff to integer powers, why does your code involve raising (sometimes-negative) stuff to non-integer powers?
That does not compute.
Regards,
Dave
MPD78
10-21-2010, 05:39 PM
Dave,
After discussing this with you and reading your last post, I see the mistake that I have made. I am going about the equation in the wrong sense, i and j are simply integer values and Cij is the matrix value.
Thanks for your help I can take it from here.
Matt
davekw7x
10-21-2010, 05:59 PM
...I can take it from here...
I can't tell you what a relief that is. I was just about running out of steam. (A little thermodynamical pun---get it? No? Oh, well...)
Regards,
Dave
MPD78
10-22-2010, 06:13 AM
I have it all working correctly now.
Thanks Dave. You are a great resource.
Matt