#include "Matrice.cpp"

#include <conio.h>
#include <stdio.h>
#include <iostream.h>
#include <stdlib.h>
#include <string.h>

typedef Matrice<double> Mat;
typedef Vecteur<double> Vect;

Mat L(3,3);

int max(int y, int x)
{
	return (x>y?x:y);
};

void Lorenz(double b, double Pr, double Ra)
{
	double temp[9]={-Pr, Pr, 0, Ra, -1, L(1,2), 0, L(2,1), -b};
        Mat tmp(3, 3, 0, temp);
	L=tmp;
};

void Lorenz(double x)
{
	L(2,3)=-x;
	L(3,2)=x;
};


void main()
{
	// Initialisations

	// Définition des paramètres:
	// b=8/3, classiquement;
	// Pr est le nombre de Prandtl, égale à 10, le plus souvent;
	// Ra=28, nombre de Rayleigh intéressant.
	int precision;
	double b,Pr,Ra,delta_t,temps;
	float b_,Pr_,Ra_,delta_t_,temps_;
	char* nom_fichier=new char[255];
	double pt0[3]={0.1,0.1,0.1};
	float pt1=0, pt2=0, pt3=0;
	Vect temp(3), k1(3), pt(3,pt0);
	double coefficients0[4]={1/4.,1/3.,1/2.,1.};
	Vect coefficients(4,coefficients0);

    printf("    Programme issu du site http://www.multimania/vmallet.\n\n");
	printf("Ce programme permet de generer des points de l'attracteur de Lorenz.\n\n");
	printf("Entrez 0 pour prendre les valeurs par defaut.\n\n");
	printf("Coefficient b (8/3 par defaut): ");
	scanf("%f",&b_);
	printf("Nombre de Prandtl (10 par defaut): ");
	scanf("%f",&Pr_);
	printf("Nombre de Rayleigh (28 par defaut): ");
	scanf("%f",&Ra_);
	printf("Condition initiale sous la forme X Y Z\n    (0 0 0 est remplace par 0.1 0.1 0.1): ");
	scanf("%f %f %f",&pt1, &pt2, &pt3);
	printf("Pas de temps (0.005 par defaut): ");
	scanf("%f",&delta_t_);
	printf("Duree de l'integration (100 par defaut): ");
	scanf("%f",&temps_);
	printf("\nPrecision des donnees (5 par defaut): ");
	scanf("%d",&precision);
	printf("Nom du fichier cible: ");
	scanf("%s",nom_fichier);

	b=b_;
	Pr=Pr_;
	Ra=Ra_;
	delta_t=delta_t_;
	temps=temps_;

	if (b==0) b=8./3.;
	if (Pr==0) Pr=10.;
	if (Ra==0) Ra=28.;
	if (delta_t==0) delta_t=0.005;
	if (temps==0) temps=100;
	if (precision==0) precision=5;
	if ((pt1!=0) || (pt2!=0) || (pt3!=0))
	{
		pt(1)=pt1;
		pt(2)=pt2;
		pt(3)=pt3;
	};

	int nb=temps/delta_t;

	printf("\nGeneration des points...\n");

	Lorenz(b,Pr,Ra);
	k1=pt;
	
	// Génération des points

	FILE *f;
	f=fopen(nom_fichier,"wb");
	char* res=new char[precision+5];
	
	for (int i=0; i<nb; i++)
	{

		// Runge-Kutta
		Lorenz(pt(1));
		for(int j=1; j<5; j++)
		{
			k1=pt+coefficients(j)*delta_t*L*k1;
		};
		pt=k1;

		_gcvt(pt(1),precision,res);
		fwrite(res,1,strlen(res),f);
		fwrite(" ",1,1,f);

		_gcvt(pt(2),precision,res);
		fwrite(res,1,strlen(res),f);
		fwrite(" ",1,1,f);

		_gcvt(pt(3),precision,res);
		fwrite(res,1,strlen(res),f);

		fwrite("\n",1,1,f);

		if (i==nb/2) printf("    Moitie des points generee...\n");
	};
}