C语言最小二乘向量旋转欧拉方法求洛伦兹方程

原创

已于 2022-01-25 09:33:46 修改 · 1.1k 阅读

5 ·

CC 4.0 BY-SA版权

文章标签：

#c语言

于 2022-01-25 09:32:50 首次发布

本文介绍了一个用于最小二乘线性拟合的C程序，包括数据读取、拟合及结果输出等功能。此外，还探讨了向量旋转及洛伦兹方程的数值解法。

最小二乘

This question will develop a set of functions to least square fit the linear model 𝑦=𝑘𝑥+𝑞 to arbitrary data provided in an input file, i.e. identify the coefficients 𝑘 and 𝑞 to optimally overlap the data points (𝑥, 𝑦) available in the input file.
本问题将开发一套函数，对输入文件中的任意数据进行线性模型𝑦=𝑘𝑥+𝑞的最小二乘拟合，即识别出系数𝑘和𝑞，使输入文件中可获得的数据点(𝑥，𝑦)最优重叠。

#include <stdio.h>
#include <stdlib.h>
#include <math.h>
typedef struct {
   
   
    char *filename;
    int num_points;
    double *x;
    double *y;
} data;
typedef struct {
   
   
double k;
double q;
double r2;
} fitcoef;
/*   4-4    */
int data_size(char filename[]){
   
   
    FILE *fp = NULL;
    int ans=0;
    double x,y;
    fp = fopen(filename, "r");
    if(fp == NULL){
   
   
        return -1;
    }
    while(fscanf(fp,"%lf %lf",&x,&y)!=EOF){
   
   
        ans++;
    }
    fclose(fp);
    return ans;
}
/*   4-5    */
int data_read(data* m_data){
   
   
    FILE *fp = NULL;
    int idx=0;
    m_data->x=(double*) malloc((m_data->num_points)*sizeof(double));
    m_data->y=(double*) malloc((m_data->num_points)*sizeof(double));
    fp = fopen(m_data->filename, "r");
    if(fp == NULL){
   
   
        return -1;
    }
    for(int i=0;i<m_data->num_points;i++){
   
   
        fscanf(fp,"%lf %lf",&(m_data->x[idx]),&(m_data->y[idx]));
        idx++;
    }
    fclose(fp);
    return 0;
}
/*   4-6    */
fitcoef linear_fit(data* m_data){
   
   
    fitcoef ans;
    double _x = 0.0, _y = 0.0;
    double s_xx = 0.0, s_yy = 0.0, s_yx = 0.0;
    for(int i=0;i<m_data->num_points;i++){
   
   
        _x+=m_data->x[i];
        _y+=m_data->y[i];
    }
    _x/=m_data->num_points;
    _y/=m_data->num_points;
    for(int i=0;i<m_data->num_points;i++){
   
   
        s_xx+=pow(m_data->x[i] - _x,2.0);
    }
    s_xx/=m_data->num_points;
    for (int i = 0; i < m_data->num_points; i++) {
   
   
        s_yy += pow((m_data->y[i] - _y),2.0) ;
    }
    s_yy/=m_data->num_points;
    for (int i = 0; i < m_data->num_points; i++) {
   
   
        s_yx += (m_data->x[i] - _x) * (m_data->y[i] - _y);
    }
    s_yx/=m_data->num_points;
    ans.k=s_yx/s_xx;
    ans.q=_y-ans.k*_x;
    ans.r2=s_yx*s_yx/s_xx/s_yy;
    return ans;
}
/*   4-7    */
void print_data(data* m_data){
   
   
    printf("\nData file: %s\n", m_data->filename);
    printf("Number of data points: %m_data\n\n", m_data->num_points   );
    printf(