Back to home page

sPhenix code displayed by LXR
 
 

    

File indexing completed on 2025-08-03 08:19:07 

0001 /*******************************************************************************
0002  * Copyright (c) 2018-2019 LongGang Pang, lgpang@qq.com
0003  *
0004  * Permission is hereby granted, free of charge, to any person obtaining a
0005  * copy of this software and/or associated documentation files (the
0006  * "Materials"), to deal in the Materials without restriction, including
0007  * without limitation the rights to use, copy, modify, merge, publish,
0008  * distribute, sublicense, and/or sell copies of the Materials, and to
0009  * permit persons to whom the Materials are furnished to do so, subject to
0010  * the following conditions:
0011  *
0012  * The above copyright notice and this permission notice shall be included
0013  * in all copies or substantial portions of the Materials.
0014  *
0015  * THE MATERIALS ARE PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND,
0016  * EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF
0017  * MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT.
0018  * IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY
0019  * CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT,
0020  * TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION WITH THE
0021  * MATERIALS OR THE USE OR OTHER DEALINGS IN THE MATERIALS.
0022  ******************************************************************************/
0023 
0024 #include <cmath>
0025 #include <cstring>
0026 #include <ctime>
0027 #include<cstdlib>
0028 #include <algorithm>
0029 #include "include/clvisc.h"
0030 #include "include/error_msgs.h"
0031 
0032 namespace clvisc {
0033 CLVisc::CLVisc(const Config & cfg, std::string device_type,
0034         int device_id):cfg_(cfg),
0035         ideal_(CLIdeal(cfg, device_type, device_id)),
0036         backend_(ideal_.get_backend()),
0037         bulkinfo_(cfg.nx, cfg.ny, cfg.nz, cfg.nxskip,
0038                 cfg.nyskip, cfg.nzskip, backend_, 
0039                 ideal_.get_compile_option())
0040 {
0041     // set current time to initial time
0042     tau0_ = cfg.tau0;
0043     tau_ = tau0_;
0044 
0045     size_ = cfg.nx * cfg.ny * cfg.nz;
0046 
0047     compile_option_ = ideal_.get_compile_option();
0048 
0049     std::cout << compile_option_ << std::endl;
0050 
0051     try {
0052         // build kernels for hydrodynamic evolution
0053         auto prg = backend_.BuildProgram("clvisc_kernel/kernel_IS.cl",
0054                 compile_option_);
0055         // is == Israel Stewart here
0056         kernel_is_initialize_ = cl::Kernel(prg, "visc_initialize");
0057         kernel_is_src_christoffel_ = cl::Kernel(prg, "visc_src_christoffel");
0058         kernel_is_src_alongx_ = cl::Kernel(prg, "visc_src_alongx");
0059         kernel_is_src_alongy_ = cl::Kernel(prg, "visc_src_alongy");
0060         kernel_is_src_alongz_ = cl::Kernel(prg, "visc_src_alongz");
0061         kernel_is_update_pimn_ = cl::Kernel(prg, "update_pimn");
0062         kernel_is_get_udiff_ = cl::Kernel(prg, "get_udiff");
0063 
0064         auto prg2 = backend_.BuildProgram("clvisc_kernel/kernel_visc.cl",
0065                 compile_option_);
0066         kernel_visc_src_christoffel_ = cl::Kernel(prg2, "kt_src_christoffel");
0067         kernel_visc_kt_src_alongx_ = cl::Kernel(prg2, "kt_src_alongx");
0068         kernel_visc_kt_src_alongy_ = cl::Kernel(prg2, "kt_src_alongy");
0069         kernel_visc_kt_src_alongz_ = cl::Kernel(prg2, "kt_src_alongz");
0070         kernel_visc_update_ev_ = cl::Kernel(prg2, "update_ev");
0071     } catch (cl::Error & err ){
0072         std::cerr<<"Error:"<<err.what()<<"("<<err.err()<<")\n";
0073         std::cerr<<"@" << __FILE__ << ":line " << __LINE__ << std::endl;
0074         std::cerr<<ErrorMessage(err.err())<<std::endl;
0075         throw(err);
0076     }
0077 }
0078 
0079 
0080 void CLVisc::initialize_gpu_buffer_() {
0081     bool read_only_option = false;
0082     // ed_bytes is the size of energy density in bytes
0083     // d_ev_[i] has size (4*ed_bytes)
0084     // d_shear_pi_[i] has size (10*ed_bytes)
0085     size_t ed_bytes = sizeof(cl_real) * size_;
0086 
0087     for (int i = 0; i < 3; i++) {
0088         d_shear_pi_[i] = backend_.CreateBufferByCopyVector(h_shear_pi_, read_only_option);
0089         d_bulk_pi_[i] = backend_.CreateBufferByCopyVector(h_bulk_pi_, read_only_option);
0090         d_net_charge_[i] = backend_.CreateBufferByCopyVector(h_net_charge_, read_only_option);
0091     }
0092 
0093     // initialize the d_udz_ to vector of (real4)
0094     // in case the nz = 1 and one wants to skip z direction calculation
0095     std::vector<cl_real4> h_udz_(size_, (cl_real4){{0.0f, 0.0f, 0.0f, 0.0f}});
0096 
0097     d_is_shear_pi_src_ = backend_.CreateBuffer(10*ed_bytes);
0098     d_is_bulk_pi_src_ = backend_.CreateBuffer(ed_bytes);
0099     d_udx_ = backend_.CreateBuffer(4*ed_bytes);    // cl_real4 type
0100     d_udy_ = backend_.CreateBuffer(4*ed_bytes);
0101     d_udz_ = backend_.CreateBufferByCopyVector(h_udz_, read_only_option);
0102     d_udiff_ = backend_.CreateBuffer(4*ed_bytes);
0103     d_goodcell_ = backend_.CreateBuffer(ed_bytes);
0104 }
0105 
0106 
0107 void CLVisc::israel_stewart_initialize_() {
0108     kernel_is_initialize_.setArg(0, d_shear_pi_[1]);
0109     kernel_is_initialize_.setArg(1, d_goodcell_);
0110     kernel_is_initialize_.setArg(2, d_udiff_);
0111     kernel_is_initialize_.setArg(3, ideal_.d_ev_[1]);
0112     kernel_is_initialize_.setArg(4, static_cast<cl_real> (tau_));
0113     kernel_is_initialize_.setArg(5, ideal_.eos_table_);
0114     backend_.enqueue_run(kernel_is_initialize_, cl::NDRange(size_), cl::NullRange);
0115 }
0116 
0117 
0118 void CLVisc::update_udiff_() {
0119     kernel_is_get_udiff_.setArg(0, d_udiff_);
0120     kernel_is_get_udiff_.setArg(1, ideal_.d_ev_[0]);
0121     kernel_is_get_udiff_.setArg(2, ideal_.d_ev_[1]);
0122     backend_.enqueue_run(kernel_is_get_udiff_,
0123             cl::NDRange(size_), cl::NullRange);
0124 }
0125 
0126 void CLVisc::half_step_israel_stewart_(int step) {
0127     kernel_is_src_christoffel_.setArg(0, d_is_shear_pi_src_);
0128     kernel_is_src_christoffel_.setArg(1, d_shear_pi_[step]);
0129     kernel_is_src_christoffel_.setArg(2, ideal_.d_ev_[step]);
0130     kernel_is_src_christoffel_.setArg(3, static_cast<cl_real> (tau_));
0131     kernel_is_src_christoffel_.setArg(4, step);
0132     backend_.enqueue_run(kernel_is_src_christoffel_,
0133             cl::NDRange(size_), cl::NullRange);
0134 
0135     kernel_is_src_alongx_.setArg(0, d_is_shear_pi_src_);
0136     kernel_is_src_alongx_.setArg(1, d_udx_);
0137     kernel_is_src_alongx_.setArg(2, d_shear_pi_[step]);
0138     kernel_is_src_alongx_.setArg(3, ideal_.d_ev_[step]);
0139     kernel_is_src_alongx_.setArg(4, ideal_.eos_table_);
0140     kernel_is_src_alongx_.setArg(5, static_cast<cl_real> (tau_));
0141     backend_.enqueue_run(kernel_is_src_alongx_, 
0142             cl::NDRange(cfg_.block_size, cfg_.ny, cfg_.nz),
0143             cl::NDRange(cfg_.block_size, 1, 1));
0144 
0145     kernel_is_src_alongy_.setArg(0, d_is_shear_pi_src_);
0146     kernel_is_src_alongy_.setArg(1, d_udy_);
0147     kernel_is_src_alongy_.setArg(2, d_shear_pi_[step]);
0148     kernel_is_src_alongy_.setArg(3, ideal_.d_ev_[step]);
0149     kernel_is_src_alongy_.setArg(4, ideal_.eos_table_);
0150     kernel_is_src_alongy_.setArg(5, static_cast<cl_real> (tau_));
0151     backend_.enqueue_run(kernel_is_src_alongy_, 
0152             cl::NDRange(cfg_.nx, cfg_.block_size, cfg_.nz),
0153             cl::NDRange(1, cfg_.block_size, 1));
0154 
0155     if (cfg_.nz != 1) {
0156         kernel_is_src_alongz_.setArg(0, d_is_shear_pi_src_);
0157         kernel_is_src_alongz_.setArg(1, d_udz_);
0158         kernel_is_src_alongz_.setArg(2, d_shear_pi_[step]);
0159         kernel_is_src_alongz_.setArg(3, ideal_.d_ev_[step]);
0160         kernel_is_src_alongz_.setArg(4, ideal_.eos_table_);
0161         kernel_is_src_alongz_.setArg(5, static_cast<cl_real> (tau_));
0162         backend_.enqueue_run(kernel_is_src_alongz_, 
0163                 cl::NDRange(cfg_.nx, cfg_.ny, cfg_.block_size),
0164                 cl::NDRange(1, 1, cfg_.block_size));
0165     }
0166 
0167     kernel_is_update_pimn_.setArg(0, d_shear_pi_[3-step]);
0168     kernel_is_update_pimn_.setArg(1, d_goodcell_);
0169     kernel_is_update_pimn_.setArg(2, d_shear_pi_[1]);
0170     kernel_is_update_pimn_.setArg(3, d_shear_pi_[step]);
0171     kernel_is_update_pimn_.setArg(4, ideal_.d_ev_[1]);
0172     kernel_is_update_pimn_.setArg(5, ideal_.d_ev_[2]);
0173     kernel_is_update_pimn_.setArg(6, d_udiff_);
0174     kernel_is_update_pimn_.setArg(7, d_udx_);
0175     kernel_is_update_pimn_.setArg(8, d_udy_);
0176     kernel_is_update_pimn_.setArg(9, d_udz_);
0177     kernel_is_update_pimn_.setArg(10, d_is_shear_pi_src_);
0178     kernel_is_update_pimn_.setArg(11, ideal_.eos_table_);
0179     kernel_is_update_pimn_.setArg(12, static_cast<cl_real> (tau_));
0180     kernel_is_update_pimn_.setArg(13, step);
0181     backend_.enqueue_run(kernel_is_update_pimn_,
0182             cl::NDRange(size_), cl::NullRange);
0183 
0184 }
0185 
0186 void CLVisc::half_step_visc_(int step) {
0187     kernel_visc_src_christoffel_.setArg(0, ideal_.d_src_);
0188     kernel_visc_src_christoffel_.setArg(1, ideal_.d_ev_[step]);
0189     kernel_visc_src_christoffel_.setArg(2, d_shear_pi_[step]);
0190     kernel_visc_src_christoffel_.setArg(3, ideal_.eos_table_);
0191     kernel_visc_src_christoffel_.setArg(4, static_cast<cl_real> (tau_));
0192     kernel_visc_src_christoffel_.setArg(5, step);
0193     backend_.enqueue_run(kernel_visc_src_christoffel_,
0194             cl::NDRange(size_), cl::NullRange);
0195 
0196     kernel_visc_kt_src_alongx_.setArg(0, ideal_.d_src_);
0197     kernel_visc_kt_src_alongx_.setArg(1, ideal_.d_ev_[step]);
0198     kernel_visc_kt_src_alongx_.setArg(2, d_shear_pi_[step]);
0199     kernel_visc_kt_src_alongx_.setArg(3, ideal_.eos_table_);
0200     kernel_visc_kt_src_alongx_.setArg(4, static_cast<cl_real> (tau_));
0201     backend_.enqueue_run(kernel_visc_kt_src_alongx_,
0202             cl::NDRange(cfg_.block_size, cfg_.ny, cfg_.nz),
0203             cl::NDRange(cfg_.block_size, 1, 1));
0204 
0205     kernel_visc_kt_src_alongy_.setArg(0, ideal_.d_src_);
0206     kernel_visc_kt_src_alongy_.setArg(1, ideal_.d_ev_[step]);
0207     kernel_visc_kt_src_alongy_.setArg(2, d_shear_pi_[step]);
0208     kernel_visc_kt_src_alongy_.setArg(3, ideal_.eos_table_);
0209     kernel_visc_kt_src_alongy_.setArg(4, static_cast<cl_real> (tau_));
0210     backend_.enqueue_run(kernel_visc_kt_src_alongy_,
0211             cl::NDRange(cfg_.nx, cfg_.block_size, cfg_.nz),
0212             cl::NDRange(1, cfg_.block_size, 1));
0213 
0214 
0215     if (cfg_.nz != 1) {
0216         kernel_visc_kt_src_alongz_.setArg(0, ideal_.d_src_);
0217         kernel_visc_kt_src_alongz_.setArg(1, ideal_.d_ev_[step]);
0218         kernel_visc_kt_src_alongz_.setArg(2, d_shear_pi_[step]);
0219         kernel_visc_kt_src_alongz_.setArg(3, ideal_.eos_table_);
0220         kernel_visc_kt_src_alongz_.setArg(4, static_cast<cl_real> (tau_));
0221         backend_.enqueue_run(kernel_visc_kt_src_alongz_,
0222                 cl::NDRange(cfg_.nx, cfg_.ny, cfg_.block_size),
0223                 cl::NDRange(1, 1, cfg_.block_size));
0224     }
0225 
0226     kernel_visc_update_ev_.setArg(0, ideal_.d_ev_[3-step]);
0227     kernel_visc_update_ev_.setArg(1, ideal_.d_ev_[1]);
0228     kernel_visc_update_ev_.setArg(2, d_shear_pi_[0]);
0229     kernel_visc_update_ev_.setArg(3, d_shear_pi_[3-step]);
0230     kernel_visc_update_ev_.setArg(4, ideal_.d_src_);
0231     kernel_visc_update_ev_.setArg(5, ideal_.eos_table_);
0232     kernel_visc_update_ev_.setArg(6, static_cast<cl_real> (tau_));
0233     kernel_visc_update_ev_.setArg(7, step);
0234     backend_.enqueue_run(kernel_visc_update_ev_,
0235             cl::NDRange(size_), cl::NullRange);
0236 } // end half_step_visc_
0237 
0238 
0239 void CLVisc::half_step_(int step) {
0240     half_step_israel_stewart_(step);
0241     half_step_visc_(step);
0242 }
0243 
0244 void CLVisc::one_step() {
0245     half_step_(1);
0246     tau_ += cfg_.dt;
0247     half_step_(2);
0248 }
0249 
0250 
0251 template <typename ValueType>
0252 void CLVisc::read_ini(const std::vector<ValueType> & ed) {
0253     ideal_.read_ini(ed);
0254     h_shear_pi_ = std::vector<cl_real> (10*size_, 0.0);
0255     h_bulk_pi_ = std::vector<cl_real> (size_, 0.0);
0256     cl_real4 zero4 = (cl_real4){{0.0, 0.0, 0.0, 0.0}};
0257     h_net_charge_ = std::vector<cl_real4> (size_, zero4);
0258     initialize_gpu_buffer_();
0259 }
0260 
0261 template <typename ValueType>
0262 void CLVisc::read_ini(const std::vector<ValueType> & ed, 
0263                       const std::vector<ValueType> & vx, 
0264                       const std::vector<ValueType> & vy, 
0265                       const std::vector<ValueType> & vz)
0266 {
0267     ideal_.read_ini(ed, vx, vy, vz);
0268     h_shear_pi_ = std::vector<cl_real> (10*size_, 0.0);
0269     h_bulk_pi_ = std::vector<cl_real> (size_, 0.0);
0270     cl_real4 zero4 = (cl_real4){{0.0, 0.0, 0.0, 0.0}};
0271     h_net_charge_ = std::vector<cl_real4> (size_, zero4);
0272     initialize_gpu_buffer_();
0273 }
0274 
0275 template <typename ValueType>
0276 void CLVisc::read_ini(const std::vector<ValueType> & ed, 
0277               const std::vector<ValueType> & vx, 
0278               const std::vector<ValueType> & vy, 
0279               const std::vector<ValueType> & vz,
0280               const std::vector<ValueType> & pi00,
0281               const std::vector<ValueType> & pi01,
0282               const std::vector<ValueType> & pi02,
0283               const std::vector<ValueType> & pi03,
0284               const std::vector<ValueType> & pi11,
0285               const std::vector<ValueType> & pi12,
0286               const std::vector<ValueType> & pi13,
0287               const std::vector<ValueType> & pi22,
0288               const std::vector<ValueType> & pi23,
0289               const std::vector<ValueType> & pi33) {
0290     ideal_.read_ini(ed, vx, vy, vz);
0291     h_shear_pi_.clear();
0292     for (size_t i = 0; i < size_; i++) {
0293         h_shear_pi_.push_back(static_cast<cl_real>(pi00.at(i)));
0294         h_shear_pi_.push_back(static_cast<cl_real>(pi01.at(i)));
0295         h_shear_pi_.push_back(static_cast<cl_real>(pi02.at(i)));
0296         h_shear_pi_.push_back(static_cast<cl_real>(pi03.at(i)));
0297         h_shear_pi_.push_back(static_cast<cl_real>(pi11.at(i)));
0298         h_shear_pi_.push_back(static_cast<cl_real>(pi12.at(i)));
0299         h_shear_pi_.push_back(static_cast<cl_real>(pi13.at(i)));
0300         h_shear_pi_.push_back(static_cast<cl_real>(pi22.at(i)));
0301         h_shear_pi_.push_back(static_cast<cl_real>(pi23.at(i)));
0302         h_shear_pi_.push_back(static_cast<cl_real>(pi33.at(i)));
0303     }
0304     h_bulk_pi_ = std::vector<cl_real> (size_, 0.0);
0305     cl_real4 zero4 = (cl_real4){{0.0, 0.0, 0.0, 0.0}};
0306     h_net_charge_ = std::vector<cl_real4> (size_, zero4);
0307     initialize_gpu_buffer_();
0308 }
0309 
0310 
0311 
0312 
0313 void CLVisc::evolve() {
0314     int max_loops = 5000;
0315     float total_exec_time = 0.0;
0316     std::time_t timer1, timer2;
0317     std::time(&timer1);
0318     tau_ = tau0_;
0319     try {
0320         israel_stewart_initialize_();
0321         ideal_.predict_first_step();
0322         for (int loop = 0; loop < max_loops; loop++) {
0323             backend_.enqueue_copy(d_shear_pi_[1], d_shear_pi_[0], 10*size_*sizeof(cl_real));
0324             backend_.enqueue_copy(ideal_.d_ev_[1], ideal_.d_ev_[0], size_*sizeof(cl_real4));
0325             one_step();
0326             update_udiff_();
0327             if (loop % cfg_.ntskip == 0) {
0328                 float max_ed = ideal_.max_energy_density();
0329                 std::cout << "tau = " << tau_ << " fm; ";
0330                 std::cout << "max_ed = " << max_ed << " ";
0331                 std::time(&timer2);
0332                 total_exec_time = std::difftime(timer2, timer1);
0333                 //std::cout << "Total computing time: " << total_exec_time << " s; ";
0334                 //std::cout << std::endl;
0335                 bulkinfo_.add_data(ideal_.d_ev_[0], d_shear_pi_[0],
0336                                    d_bulk_pi_[0], ideal_.eos_table_);
0337                 if ( max_ed < 0.05 ) break;
0338             }
0339         }
0340         std::cout << std::endl;
0341 
0342         std::cout << "Total computing time: " << total_exec_time << " s; ";
0343     } catch (cl::Error & err) {
0344         std::cout << err.what() << " " << err.err() << std::endl;
0345         std::cerr<<"@" << __FILE__ << ":line " << __LINE__ << std::endl;
0346         std::cout << ErrorMessage(err.err()) << std::endl;
0347         throw(err);
0348     }
0349 }
0350 
0351 CLVisc::~CLVisc() {
0352 }
0353 
0354 // in case the input array is a vector of double instead of float
0355 template void CLVisc::read_ini(const std::vector<float> & ed);
0356 template void CLVisc::read_ini(const std::vector<double> & ed);
0357 
0358 template void CLVisc::read_ini(const std::vector<float> & ed, 
0359                       const std::vector<float> & vx, 
0360                       const std::vector<float> & vy, 
0361                       const std::vector<float> & vz);
0362 
0363 template void CLVisc::read_ini(const std::vector<double> & ed, 
0364                       const std::vector<double> & vx, 
0365                       const std::vector<double> & vy, 
0366                       const std::vector<double> & vz);
0367 
0368 
0369 template void CLVisc::read_ini(const std::vector<float> & ed, 
0370               const std::vector<float> & vx, 
0371               const std::vector<float> & vy, 
0372               const std::vector<float> & vz,
0373               const std::vector<float> & pi00,
0374               const std::vector<float> & pi01,
0375               const std::vector<float> & pi02,
0376               const std::vector<float> & pi03,
0377               const std::vector<float> & pi11,
0378               const std::vector<float> & pi12,
0379               const std::vector<float> & pi13,
0380               const std::vector<float> & pi22,
0381               const std::vector<float> & pi23,
0382               const std::vector<float> & pi33);
0383 
0384 template void CLVisc::read_ini(const std::vector<double> & ed, 
0385               const std::vector<double> & vx, 
0386               const std::vector<double> & vy, 
0387               const std::vector<double> & vz,
0388               const std::vector<double> & pi00,
0389               const std::vector<double> & pi01,
0390               const std::vector<double> & pi02,
0391               const std::vector<double> & pi03,
0392               const std::vector<double> & pi11,
0393               const std::vector<double> & pi12,
0394               const std::vector<double> & pi13,
0395               const std::vector<double> & pi22,
0396               const std::vector<double> & pi23,
0397               const std::vector<double> & pi33);
0398 
0399 
0400 } // end namespace clvisc
[ Source navigation ] [ Diff markup ] [ Identifier search ] [ general search ]

This page was automatically generated by the 2.3.7 LXR engine.

The LXR team
Valid CSS 2.1! Valid HTML 4.01!