基于欧拉法的柯西问题近似求解C代码第5、8项功能实现求助

实现说明 #

• 第5项要求实现：首次近似计算分段数取n=1，每完成一次全区间近似计算，分段数自动加1，直到满足终止条件
• 第8项要求实现：新增小数位数计算函数，根据输入的精度eps自动确定结果需要保留的小数位数，输出时按规则舍入，同时输出迭代次数

完整修复后代码 #

#include <math.h>
#include <stdio.h>

// 按指定小数位数舍入结果
double round_y(double y, int count) {
    return round(y * pow(10, count)) / pow(10, count);
}

// 根据精度计算需要保留的小数位数（第8项核心实现）
int get_decimal_count(double eps) {
    int count = 0;
    double temp = eps;
    // 最多保留10位小数防止精度溢出
    while (temp < 1.0 && count < 10) {
        temp *= 10;
        count++;
    }
    return count;
}

// 给定的微分方程右端函数
const double f(double x, double y) {
    return 2.0 * x * sqrt(y) - y;
}

// 用指定分段数n计算区间[a,b]端点的欧拉近似解
double calc_euler(double a, double b, double s, int n, int *feasible_flag) {
    double h = (b - a) / n;
    double xi = a;
    double yi = s;
    *feasible_flag = 1;
    for (int i = 0; i < n; i++) {
        if (yi < 0) { // 可行域判断：sqrt(y)要求y>=0
            *feasible_flag = 0;
            return 0;
        }
        yi = yi + h * f(xi, yi);
        xi += h;
    }
    return yi;
}

int main() {
    int n, feasible_flag;
    double a, b, s, eps, y_prev, y_curr, error;
    printf("Enter value a = ");
    scanf("%lf", &a);
    printf("Enter value b = ");
    scanf("%lf", &b);
    printf("Enter value s (y(a)=s) = ");
    scanf("%lf", &s);
    if (s < 0) {
        printf("\nERROR Feasible region!\nPlease enter a new spacing\n");
        return 0;
    }
    printf("Enter value eps = ");
    scanf("%lf", &eps);

    // 第8项实现：计算需要保留的小数位数
    int decimal_count = get_decimal_count(eps);

    // 第5项实现：首次n=1计算
    n = 1;
    int iter_count = 0;
    y_prev = calc_euler(a, b, s, n, &feasible_flag);
    if (!feasible_flag) {
        printf("\nERROR Feasible region!\nPlease enter a new spacing\n");
        return 0;
    }

    while (1) {
        iter_count++;
        n++; // 每轮迭代分段数加1
        if (n > 10000) {
            error = fabs(y_curr - y_prev);
            printf("\nIt was not possible to achieve the specified accuracy.\n");
            printf("Current result: %.*lf\n", decimal_count, round_y(y_curr, decimal_count));
            printf("Current error: %.*lf\n", decimal_count + 2, error);
            printf("Iterations count: %d\n", iter_count);
            return 0;
        }
        y_curr = calc_euler(a, b, s, n, &feasible_flag);
        if (!feasible_flag) {
            printf("\nERROR Feasible region!\nPlease enter a new spacing\n");
            return 0;
        }
        error = fabs(y_curr - y_prev);
        if (error < eps) {
            break;
        }
        y_prev = y_curr;
    }

    // 输出结果
    printf("___________________________________\n");
    printf("Rounded final y(b) value: %.*lf\n", decimal_count, round_y(y_curr, decimal_count));
    printf("Specified accuracy eps: %.*lf\n", decimal_count + 2, eps);
    printf("Required iterations count: %d\n", iter_count);
    return 0;
}

功能验证说明 #

1. 输入边界检查：初始值s如果为负直接提示可行域错误，计算过程中出现y<0也会触发可行域提示
2. 精度控制：两次全区间近似解的差值小于eps时自动终止
3. 上限保护：n超过10000时按要求输出当前结果、误差和提示语
4. 自动舍入：结果保留的小数位数和输入的eps匹配，比如eps=0.001时结果保留3位小数

内容的提问来源于stack exchange，提问作者Anastay Kaiser