using System;
using static System.Console;
public partial class qnewton{

public static vector broyden(Func<vector,vector> f, vector x, double eps=1e-3){
	int n=x.size;
	matrix B=matrix.id(n);
	vector fx=f(x),z,fz;
	while(true){
		vector Dx=-B*fx;
x.print("broyden: x=");
fx.print("broyden: fx=");
Dx.print("broyden: Dx=");
		double s=1;
		while(true){
			z=x+Dx*s;
			fz=f(z);
			if(fz.norm()<=(1-s/2)*fx.norm()){
				Write($"good step, s={s}\n");
				break;
				}
			if(s<0.01){
				Write("bad step\n");
				B.setid();
				break;
				}
			s/=2;
		}
		vector Df=fz-fx;
		//vector c=(Dx-B*Df)/Dx.dot(Df);
		//B.update(c,Dx);
		vector c=(Dx-B*Df)/(Df%Df);
		B.update(c,Df);
B.print("B=");
		x=z;
		fx=fz;
		if(fx.norm()<eps)break;
	}
	return x;
}//broyden

}//class