numerical_integration.c

在程序设计部分主要讨论了在读写函数表达式时

#include <string.h>
#include <malloc.h>

#include "numerical_integration.h"


#ifdef _cplusplus
	extern "C" {
#endif


double sub_get_fun(double x, Token_head *head, Var_name_list *name);

void free_iter_node(Iter_node *node);


#ifdef _cplusplus
	}
#endif


int numerical_integration(char *fun, Var_name_list *name, Iter_result *result)
{

	Iter_node *data;
	Iter_node *data_head;
	
	double length;
	double prev_step;
	double cur_step;
	double cur_sum;
	double cur_correct;
	double cur_precision;
	double cur_start;
	double start = result->integration_start;
	
	double prev_t;
	double prev_sim;
	double prev_cote;
	double prev_rom;
	double preci_t;
	double preci_sim;
	double preci_cote;
	double preci_rom;


	unsigned long int iter_max = result->iter_max;
	unsigned long int iter_t = 0;
	unsigned long int iter_sim = 0;
	unsigned long int iter_cote = 0;
	unsigned long int iter_rom = 0;

	unsigned long int i;

	int return_value;

	Token_head head;

	memset(&head, 0, sizeof(Token_head));


	/*check call method*/
	if (result->status & TRAPEZOIDAL)
	{
		preci_t = fabs(result->final_precision);
		preci_sim = 0;
		preci_cote = 0;
		preci_rom = 0;

	}
	else if (result->status & SIMPSON)
	{
		preci_t = 0;
		preci_sim = fabs(result->final_precision);
		preci_cote = 0;
		preci_rom = 0;
	}
	else if (result->status & COTES)
	{
		preci_t = 0;
		preci_sim = 0;
		preci_cote = fabs(result->final_precision);
		preci_rom = 0;
	}
	else if (result->status & ROMBERG)
	{
		preci_t = 0;
		preci_sim = 0;
		preci_cote = 0;
		preci_rom = fabs(result->final_precision);
	}
	else
	{
		return ARGUMENT_FALSE;
	}

	return_value = initi_fun(fun, name, &head);
	if (return_value != OK)
	{
		return PARSE_WRONG;
	}

	data = (Iter_node *)malloc(sizeof(Iter_node));
	if (data == NULL)
	{
		return MEMERY_OVERFLOW;
	}
	//memset(data, 0, sizeof(Iter_node));
	data->preci_t = preci_t;
	data->preci_sim = preci_sim;
	data->preci_cote = preci_cote;
	data->preci_rom = preci_rom;
	data->next = NULL;
	data_head = data;

	length = result->integration_end - start;
	cur_step = length;

	data->t_value = (sub_get_fun(start, &head, name) +
					sub_get_fun(result->integration_end, &head, name)) *
						  length * 0.5;
	iter_t++;

	for (i = 0; i < iter_max; i++)
	{
		prev_step = cur_step;
		cur_step *= 0.5;
		cur_start = cur_step;
		cur_sum = sub_get_fun(start + cur_start, &head, name);

		do
		{
			cur_start += prev_step;
			if (fabs(cur_start) < fabs(length))
			{
				cur_sum += sub_get_fun(start + cur_start, &head, name);
			}
		}
		while (fabs(cur_start) < fabs(length));

		data->next = (Iter_node *)malloc(sizeof(Iter_node));
		if (data->next == NULL)
		{
			return MEMERY_OVERFLOW;
		}
		//memset(data->next, 0, sizeof(Iter_node));

		prev_t = data->t_value;
		prev_sim = data->sim_value;
		prev_cote = data->cote_value;
		prev_rom = data->rom_value;

		data = data->next;
		data->preci_t = preci_t;
		data->preci_sim = preci_sim;
		data->preci_cote = preci_cote;
		data->preci_rom = preci_rom;
		data->next = NULL;

		data->t_value = prev_t * 0.5 + cur_step * cur_sum;
		data->preci_t = fabs(prev_t - data->t_value);
		iter_t++;

		cur_correct = length;

		cur_correct *= 0.25;
		cur_precision = (data->t_value - prev_t) * length /
						(length - cur_correct);
		data->sim_value = cur_precision + prev_t;
		data->preci_sim = fabs(cur_precision);
		iter_sim++;

		if (i > 0)
		{
			cur_correct *= 0.25;
			cur_precision = (data->sim_value - prev_sim) * length /
							(length - cur_correct);
			data->cote_value = cur_precision + prev_sim;
			data->preci_cote = fabs(cur_precision);
			iter_cote++;
		}

		if (i > 1)
		{
			cur_correct *= 0.25;
			cur_precision = (data->cote_value - prev_cote) * length /
							(length - cur_correct);
			data->rom_value = cur_precision + prev_cote;
			data->preci_rom = fabs(cur_precision);
			iter_rom++;
		}


		if(data->preci_t < preci_t || data->preci_sim < preci_sim ||
		   data->preci_cote < preci_cote || data->preci_rom < preci_rom)
		{
			break;
		}
	}

	result->t_value = (double *)malloc(sizeof(double) * iter_t);
	result->t_precision = (double *)malloc(sizeof(double) * iter_t);
	result->sim_value = (double *)malloc(sizeof(double) * iter_sim);
	result->sim_precision = (double *)malloc(sizeof(double) * iter_sim);
	result->cote_value = (double *)malloc(sizeof(double) * iter_cote);
	result->cote_precision = (double *)malloc(sizeof(double) * iter_cote);
	result->rom_value = (double *)malloc(sizeof(double) * iter_rom);
	result->rom_precision = (double *)malloc(sizeof(double) * iter_rom);

	if (result->t_value == NULL || result->t_precision == NULL ||
		result->sim_value == NULL || result->sim_precision == NULL ||
		result->cote_value == NULL || result->cote_precision == NULL ||
		result->rom_value == NULL || result->rom_precision == NULL)
	{
		if (result->t_value == NULL)
		{
			free(result->t_value);
			result->t_value = NULL;
		}
		if (result->t_precision == NULL)
		{
			free(result->t_precision);
			result->t_precision = NULL;
		}
		if (result->sim_value == NULL)
		{
			free(result->sim_value);
			result->sim_value = NULL;
		}
		if (result->sim_precision == NULL)
		{
			free(result->sim_precision);
			result->sim_precision = NULL;
		}
		if (result->cote_value == NULL)
		{
			free(result->cote_value);
			result->cote_value = NULL;
		}
		if (result->cote_precision == NULL)
		{
			free(result->cote_precision);
			result->cote_precision = NULL;
		}
	    if (result->rom_value == NULL)
	    {
			free(result->rom_value);
			result->rom_value = NULL;
	    }
	    if (result->rom_precision == NULL)
	    {
			free(result->rom_precision);
			result->rom_precision = NULL;
	    }

	    return MEMERY_OVERFLOW;
	}
	
	
	data = data_head;
	result->t_value[0] = data->t_value;
	result->t_precision[0] = data->t_value;
	
	data = data->next;
	if (data != NULL)
	{
		result->t_value[1] = data->t_value;
		result->t_precision[1] = data->preci_t;
		result->sim_value[0] = data->sim_value;
		result->sim_precision[0] = data->preci_sim;
		
	
		data = data->next;
		if (data != NULL)
		{
			result->t_value[2] = data->t_value;
			result->t_precision[2] = data->preci_t;
			result->sim_value[1] = data->sim_value;
			result->sim_precision[1] = data->preci_sim;
			result->cote_value[0] = data->cote_value;
			result->cote_precision[0] = data->preci_cote;
			data = data->next;
		}
	}
	
	for (i = 3; i < iter_t && data != NULL; i++)
	{
		result->t_value[i] = data->t_value;
		result->t_precision[i] = data->preci_t;
		result->sim_value[i - 1] = data->sim_value;
		result->sim_precision[i - 1] = data->preci_sim;
		result->cote_value[i - 2] = data->cote_value;
		result->cote_precision[i - 2] = data->preci_cote;
		result->rom_value[i - 3] = data->rom_value;
		result->rom_precision[i - 3] = data->preci_rom;
		data = data->next;
	}

	free_iter_node(data_head);
	data_head = NULL;
	data = NULL;
	
	if (iter_t > 3)
	{
		result->result = result->rom_value[iter_rom - 1];
	}
	else if (iter_t > 2)
	{
		result->result = result->cote_value[iter_cote - 1];
	}
	else if (iter_t > 1)
	{
		result->result = result->sim_value[iter_sim - 1];
	}
	else
	{
		result->result = result->t_value[iter_t - 1];
	}
	
	result->total_iter = iter_t + iter_sim + iter_cote + iter_rom;
	result->iter_t = iter_t;
	result->iter_sim = iter_sim;
	result->iter_cote = iter_cote;
	result->iter_rom = iter_rom;

	free_token_head(&head);
	
	return OK;

}

double sub_get_fun(double x, Token_head *head, Var_name_list *name)
{
	double result;

	name->value = x;
	if (fun_value(head, name, &result) != OK)
	{
		return PARSE_WRONG;
	}

	return result;
}


void free_iter_node(Iter_node *node)
{
	if (node == NULL)
	{
		return;
	}

	if (node->next != NULL)
	{
		free_iter_node(node->next);
	}

	free(node);
}

void free_iter_result(Iter_result *result)
{
	if (result == NULL)
	{
		return;
	}

	if (result->t_value != NULL)
	{
		free(result->t_value);
		result->t_value = NULL;
	}
	if (result->t_precision != NULL)
	{
		free(result->t_precision);
		result->t_precision = NULL;
	}

	if (result->sim_value != NULL)
	{
		free(result->sim_value);
		result->sim_value = NULL;
	}
	if (result->sim_precision != NULL)
	{
		free(result->sim_precision);
		result->sim_precision = NULL;
	}

	if (result->cote_value != NULL)
	{
		free(result->cote_value);
		result->cote_value = NULL;
	}
	if (result->cote_precision != NULL)
	{
		free(result->cote_precision);
		result->cote_precision = NULL;
	}

	if (result->rom_value != NULL)
	{
		free(result->rom_value);
		result->rom_value = NULL;
	}
	if (result->rom_precision != NULL)
	{
		free(result->rom_precision);
		result->rom_precision = NULL;
	}

	memset(result, 0, sizeof(Iter_result));

}