Programm wont stop running

[alextmfk]

Hello!

I'm programming an Interpolation and bisection a la Lagrange
with a header-file, a main-file and a module file which I've linked.
Now I don't get any errors, but my compiled program is running to long and doesn't seem to stop. I didn't find any loop-mistakes and I'm unsure about calling a function from another function. I would very much appreciate it if anyone could point me in the right direction.

here the header:

#include<iostream>
#include<vector>
#include<fstream>
#include<cmath>

using namespace std;

double Basepol(int i,double tau, vector<double>t);

vector<double> Stuetz(int n,double t1, double t2);

double valu(double tau,vector<double>t,vector<double>z);

vector<double> Nulst(double eps,double t0,vector<double>t,vector<double>z,vector<double>x);

the not-compiled module-file:

#include "header.hpp"

double Basepol(int i,double tau,vector<double>t)
{
  double L=1.;
  for(int k=0;k<t.size();k++)
	{
		if(k!=i)
		{
			L=L*double(tau-t[k])/double(t[i]-t[k]);
		}
	}	
	return L;
}


vector<double> Stuetz(int n,double t1, double t2)
{double h=(t2-t1)/n;
 vector<double>tau;
 for (int i=0;i<=n;i++)
   {
     tau.push_back(t1+i*h);
   }
 return tau;
}

double valu(double tau,vector<double>t,vector<double>z)
{	
		double zeta=0;
		for(int i=0;i<z.size();i++)
    {	
    	zeta+=z[i]*Basepol(i,tau,t);
		}
		return zeta;
}

vector<double> Nulst(double eps,double t0,vector<double>t,vector<double>z,vector<double>x)
{
	double t1=0,tm=0,te=0;
	vector<double>erg(3,0);
	while(1)
	{	
		int j=1;
		double tmp;
		tmp=valu(t0+j,t,z);
		if(tmp<0)
		{
			t1=t0+j;
			break;
		}
		j++;
	}

	tm=(t0+t1)/2;
	while(1)
	{
	double q;
	q=valu(tm,t,z);
	if((t1-t0)<eps)break;
	if(q<0)t1=tm;
	if(q>0)t0=tm;
	tm=(t0+t1)/2;
	}
	te=(t0+t1)/2;

	erg[0]=te; erg[1]=valu(te,t,z); erg[2]=valu(te,t,x); 
	return erg;
}

and my main-file:

#include "header.hpp"

int main()
{
  vector<double>t,z,x,erg(3,0);
  double ti,zi,xi;
	int n=10;

//--------------------------------------
  ifstream istr("vertikal.dat");		
  while (istr>>ti>>zi>>xi)
	{
  	t.push_back(ti);
    z.push_back(zi);
    x.push_back(xi);
  }
  istr.close();

//--------------------------------------  
	vector<double>tau;	
	for(int i=0; i<(t.size()-1); i++)
	{
		vector<double> taui(n+1,0);
		taui=Stuetz(n,t[i],t[i+1]);
		tau.insert(tau.end(),taui.begin(),taui.end());
	}
//-------------------------------------- 
 int u=tau.size();

 	vector<double>zeta(u,0),ksi(u,0);
	
	for(int j=0;j<tau.size();j++)
  {
		for(int i=0;i<z.size();i++)
    {	
    	zeta[j]+=z[i]*Basepol(i,tau[j],t);
			ksi[j]+=x[i]*Basepol(i,tau[j],t);
    }
  }

//--------------------------------------
  ofstream ostr("Polynominterpol.dat");
  for(int i=0;i<tau.size();i++)
  	{
			ostr<<tau[i]<<"\t"<<zeta[i]<<"\t"<<ksi[i]<<endl;
    }
  ostr.close();
//--------------------------------------
	erg=Nulst(1.e-10,(t.size()-1),t,z,x);

	ofstream res("endpunkt.dat");
  res<<erg[0]<<"\t"<<erg[1]<<"\t"<<erg[2]<<endl;	
	res.close();


  return 0;
}

Thank you very much!

Alex

[lastchance]

What happens if you change

if(q<0)t1=tm;
if(q>0)t0=tm;

to

if(q<0) t1=tm;
else t0=tm;

It also would be useful debugging to cout the values of t0, t1 and q whilst this loop runs. Also, I see no evidence that you are checking that the bisection interval spans a root (function reverses sign between the two ends).

Just a piece of advice: separate the tasks of repeated bisection from the completely different task of Lagrange interpolation. Otherwise, people have to look at a wall of code, much of which is not relevant to the problem, is contained in different files (so they have to download rather than running in c++ shell), and have no access to your input files to test it.

The phrase "minimum compileable sample which illustrates the problem" sums up what you should post.

[Peter87]

I don't have the data files so I don't know where it hangs but this loop looks suspicious.

while(1)
{	
	int j=1;
	double tmp;
	tmp=valu(t0+j,t,z);
	if(tmp<0)
	{
		t1=t0+j;
		break;
	}
	j++;
}

You always pass the same values into the valu function so tmp will always get the same value. If that value is not negative there is no way it could reach the break statement.


You also need to be careful when mixing unsigned and signed integers. When I tried to run the program it crashed inside this loop because t is empty (it failed to load anything from the data files because I don't have them).

for(int i=0; i<(t.size()-1); i++)
{
	vector<double> taui(n+1,0);
	taui=Stuetz(n,t[i],t[i+1]);
	tau.insert(tau.end(),taui.begin(),taui.end());
}

t.size()-1 gives you a very large positive value because the size() function returns an unsigned integer, which means it cannot hold negative values, so it wraps around and starts counting from the largest positive value. This means the expression i<(t.size()-1) will evaluate to true when i=0 which is not so good because t[0] and t[1] are not valid indices when t is empty.

[alextmfk]

Thank you very much!

@Peter87: I defined and declared j outside the while-loop, now it works.
Of course you're right regarding signed/unsigned but it's for a course in which the professor didn't explain this yet, so I think I'll just let it be and think of you're
remark when coding in private. Thanks again for finding my bug!

@lastchance: Thanks for the advise, I'll keep it to a minimum in future :-)