I am trying to write a C++ program that finds the root of the following function:

x{(1 + [(k*n)/(1 + k*x)]} - L

using newton's method and this is the code I have and I cannot figure out what is going wrong, no matter how many iterations I put in, it never seems to converge. someone please help me!! thanks!!

#include <iostream>
#include <iomanip>
#include <math.h>
using namespace std;


enum run {TERMINATE=0, LOOP, NoRoot, HighIterations};


double f(double x, double k, double n, double l) {
return (((x*(1 + (k*n))/((1 + (k*x))))) - l);
}
double derivative(double x, double k, double n) {
return 1 + ((k*n)/((1+x)*(1+x)));
}


run newton(int NIterations, double TOL, double& x, double k, double n, double l, int summation) {


for(int iter = 1; iter <= NIterations; iter++) {
double derivativex = 1 + ((k*n)/((1+x)*(1+x)));


if(derivativex == 0.0) {
return HighIterations;
}


double deriv = -f(x, k, n, l)/derivativex;
x = x + deriv;


if(fabs(deriv) < TOL) {
return TERMINATE;
}


cout << "Iteration " << iter << ": x = " << x << ", dx = " << deriv << endl;
}
return LOOP;


double sum = 0;
for (int i = 0; i < summation + 1; i++) {   //perform the function
sum = sum + (x*(1 + (k*n)/(1 + (k*x))) - l);
}
}


int main() {
double zero;
double tolerance;
double k;
double n;
double l;
int iterations;
int summ;


cout << "please enter a value for M, the number of different binding molecules: ";
cin >> summ;
for (int q = 0; q < summ; q++) {
cout << "please enter a value for k, the equilibrium constant(s): ";
cin >> k;
}
for (int z = 0; z < summ; z++) {
cout << "please enter a value for n, the total concentration(s) of the various binding molecules: ";
cin >> n;
}
cout << "please enter value, l, for the total concentration of the ligand: ";
cin >> l;
cout << "please enter a value that you may guess to be the root: ";
cin >> zero;
cout << "please enter a value for the tolerance level: ";
cin >> tolerance;
if (tolerance <= 0) {
cout << "The tolerance level must be greater than zero.  Please enter a new value: ";
cin >> tolerance;
}
cout << "please enter the numer of iterations you would like to perform: ";
cin >> iterations;


run w = newton(iterations, tolerance, zero, k, n, l, summ);


switch(w) {
case TERMINATE: {
cout << "The root is " << zero << endl;
return 0;
}
case LOOP:
cout << "FAIL: Did not converge in " << iterations << " iterations!" << endl;
}
}
Member Avatar for iamthwee

I've not tested your code, but there could be a variety of problems there.

1. Are you sure that the derivative of the function

x*(1+(k*n/(1+k*x)))-l

too lazy to Latex.it.up

is correct?

Me thinks not. Check this out...

http://www.calc101.com/webMathematica/derivatives.jsp#topdoit

And on entering the function exactly as shown below:-

x*(1+(k*n/(1+k*x)))-l

yields the derivative as being:-

[(k^2)(x^2) + k(n+2x) +1 ]/ (kx+1)^2

-In fact the equation you have given is not a trivial one to solve. Application of the chain rule and other crapola is needed.

Second, is that you function has a complex shape. You need to ensure you choose the correct start value otherwise it may not converge.

Anyway, my guess would be that your derivative is wrong. So check that out first! Ho ho ho.

use the

[ Code ]
and
[ /Code ]


tags to post ur code...it helps out alot and also take the spaces out of the [C like that

Use CODE tags

Better late than never:

#include <iostream>
#include <iomanip>
#include <math.h>
using namespace std;

enum run {TERMINATE=0, LOOP, NoRoot, HighIterations};

double f(double x, double k, double n, double l) {
return (((x*(1 + (k*n))/((1 + (k*x))))) - l);
}
double derivative(double x, double k, double n) { 
return 1 + ((k*n)/((1+x)*(1+x)));
}

run newton(int NIterations, double TOL, double& x, double k, double n, double l, int summation) {

for(int iter = 1; iter <= NIterations; iter++) {
double derivativex = 1 + ((k*n)/((1+x)*(1+x)));

if(derivativex == 0.0) {
return HighIterations;
}

double deriv = -f(x, k, n, l)/derivativex;
x = x + deriv;

if(fabs(deriv) < TOL) {
return TERMINATE;
}

cout << "Iteration " << iter << ": x = " << x << ", dx = " << deriv << endl;
} 
return LOOP;

double sum = 0;
for (int i = 0; i < summation + 1; i++) { //perform the function
sum = sum + (x*(1 + (k*n)/(1 + (k*x))) - l);
}
}

int main() {
double zero;
double tolerance;
double k;
double n;
double l;
int iterations;
int summ;

cout << "please enter a value for M, the number of different binding molecules: ";
cin >> summ;
for (int q = 0; q < summ; q++) {
cout << "please enter a value for k, the equilibrium constant(s): ";
cin >> k;
}
for (int z = 0; z < summ; z++) {
cout << "please enter a value for n, the total concentration(s) of the various binding molecules: ";
cin >> n;
}
cout << "please enter value, l, for the total concentration of the ligand: ";
cin >> l;
cout << "please enter a value that you may guess to be the root: ";
cin >> zero;
cout << "please enter a value for the tolerance level: ";
cin >> tolerance;
if (tolerance <= 0) {
cout << "The tolerance level must be greater than zero. Please enter a new value: ";
cin >> tolerance;
}
cout << "please enter the numer of iterations you would like to perform: ";
cin >> iterations;

run w = newton(iterations, tolerance, zero, k, n, l, summ);

switch(w) {
case TERMINATE: {
cout << "The root is " << zero << endl;
return 0;
}
case LOOP:
cout << "FAIL: Did not converge in " << iterations << " iterations!" << endl;
}
}
Be a part of the DaniWeb community

We're a friendly, industry-focused community of developers, IT pros, digital marketers, and technology enthusiasts meeting, networking, learning, and sharing knowledge.