omp/solve.cpp - master - NenNil/vmk2025-parallel-computing

omp/solve.cpp
167 строк5 KB
Нил Ненахов
optimization
02 ноя 2025, 00:22
02 ноя 2025, 00:2215805f8
#include "solve.h"
#include "matrix.h"
#include <cmath>
#include <iostream>
#include <mpi.h>
#include <string>
#include <vector>
#include <chrono>

using namespace std;

void printVector(const std::vector<double>& vec) {
    for (const auto& element : vec) { // Используем цикл range-based for
        std::cout << element << ' ';
    }
    std::cout << '\n';
}

// Решение системы линейных уравнений методом сопряженных градиентов
void solve(CSRMatrix* matrix, const std::vector<double>& b, double eps, int maxit,
           std::vector<double>& x, int* n, std::vector<double>& z,
           std::vector<double>& p, std::vector<double>& q, CSRMatrix* M) {

    int nc = matrix->rows;
    // Инициализация вектора x нулями
    // vfill(nc, x, 0.0);

    std::vector<double> r = b;
    double rho = 0.0, prho = 0.0, beta = 0.0, alpha = 0.0;
    int k = 0;

    std::vector<double> h(r.size());
    double check;

    // Копирование диагонали из matrix в M
    // copy_diag(M, matrix);

    std::chrono::time_point<std::chrono::high_resolution_clock> start, end;
    double duration;   
    start = std::chrono::high_resolution_clock::now();


    do {

    end = std::chrono::high_resolution_clock::now();
    duration = std::chrono::duration_cast<std::chrono::microseconds>(end - start).count();
    duration = duration/1000000;
    std::cout << "101:" << duration << ";";
    

    start = std::chrono::high_resolution_clock::now();

        ++k;

        // Вычисление SpMV: z = M * r
        SpMV(M, r, z);


    end = std::chrono::high_resolution_clock::now();
    duration = std::chrono::duration_cast<std::chrono::microseconds>(end - start).count();
    duration = duration/1000000;
    std::cout << "102:" << duration << ";";
    

    start = std::chrono::high_resolution_clock::now();

        // Обновление rho: rho = dot(r, z)
        prho = rho;
        rho = dot(nc, r, z);


    end = std::chrono::high_resolution_clock::now();
    duration = std::chrono::duration_cast<std::chrono::microseconds>(end - start).count();
    duration = duration/1000000;
    std::cout << "103:" << duration << ";";
    

    start = std::chrono::high_resolution_clock::now();

        if (k == 1) {
            // Первый шаг: p = z
            vcopy(nc, z, p);
        } else {
            // Вычисление beta: beta = rho / prho
            beta = rho / prho;

            // Обновление p: p = z + beta * p
            axpy(nc, z, p, beta, p);
        }


    end = std::chrono::high_resolution_clock::now();
    duration = std::chrono::duration_cast<std::chrono::microseconds>(end - start).count();
    duration = duration/1000000;
    std::cout << "104:" << duration << ";";
    

    start = std::chrono::high_resolution_clock::now();

        // Вычисление SpMV: q = matrix * p
        SpMV(matrix, p, q);


    end = std::chrono::high_resolution_clock::now();
    duration = std::chrono::duration_cast<std::chrono::microseconds>(end - start).count();
    duration = duration/1000000;
    std::cout << "105:" << duration << ";";
    

    start = std::chrono::high_resolution_clock::now();

        // Вычисление alpha: alpha = rho / dot(p, q)
        alpha = rho / dot(matrix->IA.size()-1, p, q);


    end = std::chrono::high_resolution_clock::now();
    duration = std::chrono::duration_cast<std::chrono::microseconds>(end - start).count();
    duration = duration/1000000;
    std::cout << "106:" << duration << ";";
    

    start = std::chrono::high_resolution_clock::now();

        // Обновление x: x = x + alpha * p
        axpy(nc, x, p, alpha, x);


    end = std::chrono::high_resolution_clock::now();
    duration = std::chrono::duration_cast<std::chrono::microseconds>(end - start).count();
    duration = duration/1000000;
    std::cout << "107:" << duration << ";";
    

    start = std::chrono::high_resolution_clock::now();

        // Обновление r: r = r - alpha * q
        axpy(nc, r, q, -alpha, r);


    end = std::chrono::high_resolution_clock::now();
    duration = std::chrono::duration_cast<std::chrono::microseconds>(end - start).count();
    duration = duration/1000000;
    std::cout << "108:" << duration << ";";
    

    start = std::chrono::high_resolution_clock::now();

        /*
        //Проверка
        cout << "J(w) = ";
        for (int i = 0; i<r.size(); i++)
        {
            h[i] = r[i] + b[i]; 
        }
        cout << " " << -0.5*dot(r.size(), h, x) << endl;
        */

    } while ((rho > eps * eps) && (k < maxit));

    *n = k;

    end = std::chrono::high_resolution_clock::now();
    duration = std::chrono::duration_cast<std::chrono::microseconds>(end - start).count();
    duration = duration/1000000;
    std::cout << "109:" << duration << ";";
    
}