/*
 * Pseudospectral 2nd order 2D linear PDE solver
 * Copyright (C) 2016 Anton Khirnov <anton@khirnov.net>
 *
 * This program is free software: you can redistribute it and/or modify
 * it under the terms of the GNU General Public License as published by
 * the Free Software Foundation, either version 3 of the License, or
 * (at your option) any later version.
 *
 * This program is distributed in the hope that it will be useful,
 * but WITHOUT ANY WARRANTY; without even the implied warranty of
 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
 * GNU General Public License for more details.
 *
 * You should have received a copy of the GNU General Public License
 * along with this program.  If not, see <http://www.gnu.org/licenses/>.
 */

#ifndef QMS_PSSOLVE_H
#define QMS_PSSOLVE_H

#include "common.h"

#if HAVE_OPENCL
#include <cl.h>
#else
typedef void* cl_context;
typedef void* cl_command_queue;
#endif

#include <stdint.h>

#include "basis.h"

enum PSSolveDiffOrder {
    PSSOLVE_DIFF_ORDER_00,
    PSSOLVE_DIFF_ORDER_10,
    PSSOLVE_DIFF_ORDER_01,
    PSSOLVE_DIFF_ORDER_11,
    PSSOLVE_DIFF_ORDER_20,
    PSSOLVE_DIFF_ORDER_02,
    PSSOLVE_DIFF_ORDER_NB,
};

typedef struct PSSolvePriv PSSolvePriv;

typedef struct PSSolveContext {
    /**
     * Solver private data, not to be touched by the caller.
     */
    PSSolvePriv *priv;

    /**
     * The basis sets to be used in each direction.
     * Set by the caller before qms_pssolve_context_init().
     */
    const BasisSet *basis[2];

    /**
     * Order of the solver in each direction.
     * Set by the caller before qms_pssolve_context_init().
     */
    int solve_order[2];

    /**
     * Locations of the collocation points in each direction. The equation
     * coefficients passed to qms_pssolve_solve() should be evaluated at those
     * grid positions.
     *
     * Set by the solver after qms_pssolve_context_init(). Each array is
     * solve_order[i]-sized.
     */
    double *colloc_grid[2];

    cl_context       ocl_ctx;
    cl_command_queue ocl_queue;

    uint64_t lu_solves_count;
    uint64_t lu_solves_time;

    uint64_t cg_solve_count;
    uint64_t cg_iter_count;
    uint64_t cg_time_total;

    uint64_t construct_matrix_count;
    uint64_t construct_matrix_time;
} PSSolveContext;

/**
 * Allocate a new solver.
 */
int qms_pssolve_context_alloc(PSSolveContext **ctx);

/**
 * Initialize the solver for use after all the context options have been set.
 */
int qms_pssolve_context_init(PSSolveContext *ctx);

/**
 * Free the solver and all its internal state.
 */
void qms_pssolve_context_free(PSSolveContext **ctx);

/**
 * Solve a second order linear PDE in 2D with a pseudospectral method.
 *
 * @param eq_coeffs the coefficients of each derivative term at the collocation
 *                  points.
 * @param rhs the right-hand side of the equation at the collocation points.
 * @param coeffs the spectral coefficients of the solution will be written here.
 */
int qms_pssolve_solve(PSSolveContext *ctx,
                      const double * const eq_coeffs[PSSOLVE_DIFF_ORDER_NB],
                      const double *rhs, double *coeffs);

#endif /* QMS_PSSOLVE_H */