Added solvers

This commit is contained in:
Alberto Ramos 2022-05-31 13:10:47 +02:00
parent 6d965505ec
commit a5e7cfd8b0
2 changed files with 83 additions and 0 deletions

55
src/Solvers/CG.jl Normal file
View file

@ -0,0 +1,55 @@
###
### "THE BEER-WARE LICENSE":
### Alberto Ramos wrote this file. As long as you retain this
### notice you can do whatever you want with this stuff. If we meet some
### day, and you think this stuff is worth it, you can buy me a beer in
### return. <alberto.ramos@cern.ch>
###
### file: CG.jl
### created: Tue Nov 30 11:10:57 2021
###
"""
function CG!
Solves the linear equation `Ax = si`
"""
function CG!(si, U, m0, A, lp::SpaceParm, dws::DiracWorkspace)
dws.sr .= si
dws.sp .= si
norm = CUDA.mapreduce(x -> norm2(x), +, si)
fill!(si,zero(eltype(so)))
err = 0.0
tol = eps * norm
iterations = 0
niter = 0
for i in 1:maxiter
A(dws.sAp, U, dws.sp, am0, dws.st, lp)
prod = CUDA.mapreduce(x -> dot(x[1],x[2]), +, zip(dws.sp, dws.sAp))
alpha = norm/prod
si .= si .+ alpha .* dws.sp
dws.sr .= dws.sr .- alpha .* dws.sAp
err = CUDA.mapreduce(x -> norm2(x), +, dws.sr)
if err < tol
niter = i
break
end
beta = err/norm
dws.sp .= dws.sr .+ beta .* dws.sp
norm = err;
end
if err > tol
error("CG! not converged after $maxiter iterations (Residuals: $err)")
end
return niter
end

28
src/Solvers/Solvers.jl Normal file
View file

@ -0,0 +1,28 @@
###
### "THE BEER-WARE LICENSE":
### Alberto Ramos wrote this file. As long as you retain this
### notice you can do whatever you want with this stuff. If we meet some
### day, and you think this stuff is worth it, you can buy me a beer in
### return. <alberto.ramos@cern.ch>
###
### file: Solvers.jl
### created: Tue Nov 30 11:09:51 2021
###
module Solvers
using CUDA, TimerOutputs
using ..Space
using ..Groups
using ..Fields
using ..YM
using ..Spinors
using ..Dirac
include("CG.jl")
export CG!
end