anlt_wind_solver_opencl.cc

#include <vtestbed/config/openmp.hh>
#include <vtestbed/config/opencl.hh>
#include <vtestbed/config/real_type.hh>
#include <vtestbed/core/anlt_wind_solver.hh>
#include <vtestbed/linalg/linear_algebra.hh>
#include <vtestbed/core/ship_hull_panel.hh>
#include <vtestbed/core/policy.hh>
#include <vtestbed/opencl/opencl.hh>
#include <vtestbed/opencl/vector.hh>

using namespace vtb::opencl;

namespace vtb {

    namespace core {

        template <class T>
        class ANLT_wind_solver_opencl:
            public ANLT_wind_solver<T>,
            public Context_base {

        private:
            using base_type = ANLT_wind_solver<T>;
            using typename base_type::array_type;
            using typename base_type::wind_field_type;
            using typename base_type::grid_type;
            using typename base_type::panel_type;
            using typename base_type::panel_array;
            using typename base_type::ray_type;
            using typename base_type::ray_array;
            using vec3 = blitz::TinyVector<T,3>;

        private:
            clx::kernel _kernel_near_the_boundary, _kernel_on_the_boundary;
            Buffer<vec3> _d_wind_field;

        public:

            using Context_base::context;

            void context(Context* rhs) override {
                Context_base::context(rhs);
                auto prog = context()->compiler().compile("anlt_wind.cl");
                this->_kernel_near_the_boundary = prog.kernel("wind_near_the_boundary");
                this->_kernel_on_the_boundary = prog.kernel("wind_on_the_boundary");
            }

            void step(const grid_type& grid, panel_array& panels) override;

        };
    }
}

template <class T>
void vtb::core::ANLT_wind_solver_opencl<T>::step(
    const grid_type& grid,
    panel_array& panels
) {
    const int npanels = panels.size();
    auto& wind_field = this->_wind_field;
    if (npanels == 0) {
        if (this->near_the_boundary()) {
            wind_field.resize(grid.shape());
            wind_field = T{};
        }
        return;
    }
    auto& ppl = this->context()->pipeline();
    vec3 u(5,0,0);
    Buffer<panel_type> d_panels;
    ppl.copy(panels, d_panels);
    if (this->near_the_boundary()) {
        auto& d_wind_field = this->_d_wind_field;
        wind_field.resize(grid.shape());
        ppl.allocate(wind_field, d_wind_field);
        auto& kernel = this->_kernel_near_the_boundary;
        kernel.argument(0, npanels);
        kernel.argument(1, d_panels);
        kernel.argument(2, d_wind_field);
        kernel.argument(3, make_vector(grid.lbound()));
        kernel.argument(4, make_vector(grid.ubound()));
        kernel.argument(5, make_vector(grid.shape()));
        kernel.argument(6, make_vector(u));
        ppl.step();
        ppl.kernel(kernel, clx::range(wind_field.size()));
        ppl.step();
        ppl.copy(d_wind_field, wind_field);
    }
    if (this->on_the_boundary()) {
        auto& kernel = this->_kernel_on_the_boundary;
        kernel.argument(0, npanels);
        kernel.argument(1, d_panels);
        kernel.argument(2, make_vector(u));
        ppl.step();
        ppl.kernel(kernel, clx::range(npanels));
        ppl.step();
        ppl.copy(d_panels, panels);
    }
    if (this->on_the_boundary() || this->near_the_boundary()) {
        ppl.wait();
    }
}

template <>
auto
vtb::core::make_anlt_wind_solver<VTB_REAL_TYPE,vtb::core::Policy::OpenCL>() ->
std::unique_ptr<ANLT_wind_solver<VTB_REAL_TYPE>> {
    return std::unique_ptr<ANLT_wind_solver<VTB_REAL_TYPE>>(
        new ANLT_wind_solver_opencl<VTB_REAL_TYPE>
    );
}

template class vtb::core::ANLT_wind_solver_opencl<VTB_REAL_TYPE>;