From 9a21acd651338c0ed5602db448d7a52e9fcc98d3 Mon Sep 17 00:00:00 2001 From: houpc Date: Thu, 11 Apr 2024 11:24:01 +0800 Subject: [PATCH 1/3] update diagsGV() using Graph AND remove framework argument in python compiler --- README.md | 4 +- src/backend/compiler_python.jl | 20 +-- src/frontend/GV.jl | 20 ++- src/frontend/GV_diagrams/readfile.jl | 195 +++++++++++++++++++++++++-- 4 files changed, 209 insertions(+), 30 deletions(-) diff --git a/README.md b/README.md index 4d9e3540..c3063ad2 100644 --- a/README.md +++ b/README.md @@ -111,8 +111,8 @@ Compilers.compile_Julia(g_o21, "func_o21.jl"); # Export the self-energy's source code to a C file. Compilers.compile_C(g_o21, "func_o21.c"); -# Export the self-energy's source code to a Python file for use with the JAX machine learning framework. -Compilers.compile_Python(g_o21, :jax, "func_o21.py"); +# Export the self-energy's source code to a Python file. +Compilers.compile_Python(g_o21, "func_o21.py"); ``` ### Computational Graph visualization diff --git a/src/backend/compiler_python.jl b/src/backend/compiler_python.jl index 3fd231e1..1efc40c9 100644 --- a/src/backend/compiler_python.jl +++ b/src/backend/compiler_python.jl @@ -1,20 +1,13 @@ """ function to_python_str(graphs::AbstractVector{<:AbstractGraph}) - Compile a list of graphs into a string for a python static function and output a python script which support the mindspore and jax framework. + Compile a list of graphs into a string for a python static function and output a python script. # Arguments: - `graphs` vector of computational graphs -- `framework` the type of the python frameworks, including `:jax` and `mindspore`. """ -function to_python_str(graphs::AbstractVector{<:AbstractGraph}, framework::Symbol=:jax) - if framework == :jax - head = "from jax import jit\n" - elseif framework == :mindspore - head = "import mindspore as ms\n@ms.jit\n" - else - error("no support for $type framework") - end +function to_python_str(graphs::AbstractVector{<:AbstractGraph}) + head = "" body = "" leafidx = 0 root = [id(g) for g in graphs] @@ -52,15 +45,12 @@ function to_python_str(graphs::AbstractVector{<:AbstractGraph}, framework::Symbo output = join(output, ",") tail = " return $output\n\n" - if framework == :jax - tail *= "graphfunc_jit = jit(graphfunc)" - end expr = head * body * tail return expr, gid_to_leafid end -function compile_Python(graphs::AbstractVector{<:AbstractGraph}, framework::Symbol=:jax, filename::String="GraphFunc.py") - py_string, leafmap = to_python_str(graphs, framework) +function compile_Python(graphs::AbstractVector{<:AbstractGraph}, filename::String="GraphFunc.py") + py_string, leafmap = to_python_str(graphs) open(filename, "w") do f write(f, py_string) end diff --git a/src/frontend/GV.jl b/src/frontend/GV.jl index a9831670..762409ef 100644 --- a/src/frontend/GV.jl +++ b/src/frontend/GV.jl @@ -16,7 +16,7 @@ import ..FrontEnds: Proper #ver4, ver3, and polarization diagrams may require t import ..FrontEnds: Response, Composite, ChargeCharge, SpinSpin, UpUp, UpDown import ..FrontEnds: AnalyticProperty, Instant, Dynamic import ..FrontEnds: TwoBodyChannel, Alli, PHr, PHEr, PPr, AnyChan -import ..FrontEnds: DiagramId, Ver4Id, Ver3Id, GreenId, SigmaId, PolarId, BareGreenId, BareInteractionId +import ..FrontEnds: DiagramId, Ver4Id, Ver3Id, GreenId, SigmaId, PolarId, GenericId, BareGreenId, BareInteractionId using AbstractTrees @@ -74,6 +74,24 @@ function diagsGV(type::Symbol, order::Int, GOrder::Int=0, VerOrder::Int=0; end end +function diagsGV(type::Symbol, order::Int; spinPolarPara::Float64=0.0, filter::Vector{Filter}=[NoHartree]) + if type == :spinPolar + filename = string(@__DIR__, "/GV_diagrams/groups_spin/Polar$(order)_0_0.diag") + elseif type == :chargePolar + filename = string(@__DIR__, "/GV_diagrams/groups_charge/Polar$(order)_0_0.diag") + elseif type == :sigma + filename = string(@__DIR__, "/GV_diagrams/groups_sigma/Sigma$(order)_0_0.diag") + elseif type == :green + filename = string(@__DIR__, "/GV_diagrams/groups_green/Green$(order)_0_0.diag") + elseif type == :freeEnergy + filename = string(@__DIR__, "/GV_diagrams/groups_free_energy/FreeEnergy$(order)_0_0.diag") + else + error("no support for $type diagram") + end + + return read_diagrams(filename, type, filter=filter, spinPolarPara=spinPolarPara) +end + """ function diagsGV_ver4(order::Int; spinPolarPara::Float64=0.0, filter::Vector{Filter}=[NoHartree], channels=[PHr, PHEr, PPr, Alli]) diff --git a/src/frontend/GV_diagrams/readfile.jl b/src/frontend/GV_diagrams/readfile.jl index 01c91ae9..f7fbcef1 100644 --- a/src/frontend/GV_diagrams/readfile.jl +++ b/src/frontend/GV_diagrams/readfile.jl @@ -68,17 +68,17 @@ function getK(loopNum::Int, loopIdx::Int) return k end -function firstLoopIdx(type) +function extK_num(type) if type == :vertex4 || type == :vertex4I #three extK - return 4 + return 3 elseif type == :sigma #one extK - return 2 + return 1 elseif type == :green #one extK - return 2 + return 1 elseif type == :chargePolar || type == :spinPolar #one extK - return 2 - elseif type == :freeEnergy #no extK return 1 + elseif type == :freeEnergy #no extK + return 0 else error("not implemented!") end @@ -189,7 +189,7 @@ function read_diagrams(filename::AbstractString; labelProd::Union{Nothing,LabelP end function read_vertex4diagrams(filename::AbstractString; spinPolarPara::Float64=0.0, filter=[NoHartree], channels=[PHr, PHEr, PPr, Alli], - keywords::Vector{String}=["SelfEnergy", "DiagNum", "Order", "GNum", "Ver4Num", "LoopNum", "ExtLoopIndex", + keywords::Vector{String}=["Vertex4", "DiagNum", "Order", "GNum", "Ver4Num", "LoopNum", "ExtLoopIndex", "DummyLoopIndex", "TauNum", "DummyTauIndex"]) # Open a diagram file io = open(filename, "r") @@ -248,8 +248,6 @@ function read_vertex4diagrams(filename::AbstractString; spinPolarPara::Float64=0 unique!(gr_keys) graphvec = Graph{_dtype.factor,_dtype.weight}[] - # extTs = Vector{NTuple{4,Int}}() - # responses = Vector{Response}() for (extT, channel) in gr_keys keyDi = (extT, channel, 0) keyEx = (extT, channel, 1) @@ -261,10 +259,8 @@ function read_vertex4diagrams(filename::AbstractString; spinPolarPara::Float64=0 guuId = Ver4Id(para, UpUp, gId_Di.type, k=gId_Di.extK, t=gId_Di.extT, chan=gId_Di.channel) guu = Graph([gud, gEx], properties=guuId) append!(graphvec, [guu, gud]) - # append!(extTs, [extT, extT]) - # append!(responses, [UpUp, UpDown]) end - # return graphvec, extTs, responses + return graphvec end @@ -413,6 +409,181 @@ function read_one_vertex4diagram!(io::IO, GNum::Int, verNum::Int, loopNum::Int, return g_Di, g_Ex end +function read_diagrams(filename::AbstractString, diagType::Symbol; filter=[NoHartree], spinPolarPara::Float64=0.0, + keywords::Vector{String}=["SelfEnergy", "DiagNum", "Order", "GNum", "Ver4Num", "LoopNum", "ExtLoopIndex", + "DummyLoopIndex", "TauNum", "ExtTauIndex", "DummyTauIndex"] +) + # Open a diagram file + io = open(filename, "r") + + # Read global graph properties + diagNum, loopNum, tauNum, verNum = 1, 1, 2, 0 + extIndex = Int[] + GNum = 2 + lineNum = 1 + while true + line = readline(io) + length(line) == 0 && break + keyword = keywords[lineNum] + # @assert occursin(keyword, line) + if keyword == "DiagNum" + diagNum = _StringtoIntVector(line)[1] + elseif keyword == "GNum" + GNum = _StringtoIntVector(line)[1] + elseif keyword == "Ver4Num" + verNum = _StringtoIntVector(line)[2] + elseif keyword == "LoopNum" + loopNum = _StringtoIntVector(line)[1] + elseif keyword == "TauNum" + tauNum = _StringtoIntVector(line)[1] + elseif keyword == "ExtTauIndex" + extIndex = _StringtoIntVector(line) + end + lineNum += 1 + end + offset_ver4 = diagType == :sigma ? 1 : 0 + + # Read one diagram at a time + diagrams = Graph{_dtype.factor,_dtype.weight}[] + for _ in 1:diagNum + diag = read_one_diagram!(diagType, IOBuffer(readuntil(io, "\n\n")), + GNum, verNum, loopNum, extIndex, spinPolarPara, filter=filter, offset_ver4=offset_ver4) + # isempty(diags) && continue + push!(diagrams, diag) + end + close(io) + + if diagType == :freeEnergy + return [IR.linear_combination(diagrams, properties=diagrams[1].properties),] + else + extT_labels = Vector{NTuple{length(extIndex),Int}}() + sizehint!(extT_labels, length(diagrams)) + for prop in IR.properties.(diagrams) + push!(extT_labels, prop.extT) + end + gr = _group(diagrams, extT_labels) + unique!(extT_labels) + graphvec = Graph{_dtype.factor,_dtype.weight}[] + for key in extT_labels + gId = gr[key][1].properties + push!(graphvec, IR.linear_combination(gr[key], properties=gId)) + end + return graphvec + end +end + +function read_one_diagram!(type::Symbol, io::IO, GNum::Int, verNum::Int, loopNum::Int, extIndex::Vector{Int}, + spinPolarPara::Float64=0.0; filter=[NoHartree], splitter="|", offset::Int=-1, offset_ver4::Int=0 +) + flag_proper = false + if Proper in filter + flag_proper = true + end + isDynamic = verNum == 1 ? false : true + ################ Read Hugenholtz Diagram information #################### + @assert occursin("Permutation", readline(io)) + permutation = _StringtoIntVector(readline(io)) .- offset + @assert length(permutation) == length(unique(permutation)) == GNum + + @assert occursin("SymFactor", readline(io)) + symfactor = parse(Float64, readline(io)) + + @assert occursin("GType", readline(io)) + opGType = _StringtoIntVector(readline(io)) + @assert length(opGType) == GNum + + @assert occursin("VertexBasis", readline(io)) + tau_labels = _StringtoIntVector(readline(io)) .- offset + readline(io) + + @assert occursin("LoopBasis", readline(io)) + currentBasis = zeros(Int, (GNum, loopNum)) + for i in 1:loopNum + x = parse.(Int, split(readline(io))) + @assert length(x) == GNum + currentBasis[:, i] = x + end + + @assert occursin("Ver4Legs", readline(io)) + if verNum == 0 + ver4Legs = Vector{Vector{Int64}}(undef, 0) + else + strs = split(readline(io), splitter) + ver4Legs = _StringtoIntVector.(strs[1:verNum]) + end + + @assert occursin("WType", readline(io)) + if verNum > 0 + opWType = _StringtoIntVector(readline(io)) + end + + @assert occursin("SpinFactor", readline(io)) + spinFactors = _StringtoIntVector(readline(io)) + + extIndex = extIndex .- offset + extNum = length(extIndex) + extK = zeros(loopNum) + + greens = Graph{_dtype.factor,_dtype.weight}[] + # create all fermionic propagators + for (ind1, ind2) in enumerate(permutation) + opGType[ind1] == -2 && continue + diagid = BareGreenId(k=currentBasis[ind1, :], t=[tau_labels[ind1], tau_labels[ind2]]) + push!(greens, Graph([]; properties=diagid)) + end + fermi_greenProd = Graph(greens, operator=IR.Prod()) + + interactions = Graph{_dtype.factor,_dtype.weight}[] + spinfactors_existed = Float64[] + # println("##### $permutation $ver4Legs") + for (iex, spinFactor) in enumerate(spinFactors) + # create permutation and ver4Legs for each Feynman diagram from a Hugenholtz diagram + spinFactor == 0 && continue + # flag_proper && proper[iex] == 1 && continue + push!(spinfactors_existed, sign(spinFactor) * (2 / (1 + spinPolarPara))^(log2(abs(spinFactor)))) + + permu, ver4Legs_ex = _exchange(permutation, ver4Legs, iex, extNum, offset_ver4=offset_ver4) + + ######################## Create Feynman diagram ######################### + leafs = Graph{_dtype.factor,_dtype.weight}[] + + # create all bosionic operators (relevant to interaction lines) + for verLeg in ver4Legs_ex + ind1, ind2 = verLeg[2] - offset, verLeg[4] - offset + current = currentBasis[verLeg[1]-offset, :] - currentBasis[ind1, :] + @assert current == currentBasis[ind2, :] - currentBasis[verLeg[3]-offset, :] # momentum conservation + + diagid = BareInteractionId(ChargeCharge, k=current, t=[tau_labels[ind1], tau_labels[ind2]]) + push!(leafs, Graph([]; properties=diagid)) + end + isempty(leafs) && continue + # push!(interactions, Graph(leafs, operator=IR.Prod(), factor=spinfactor * symfactor)) + push!(interactions, Graph(leafs, operator=IR.Prod())) + end + + innerLoopNum = loopNum - extNum + 1 + if type == :freeEnergy + innerLoopNum -= 1 + diagid = GenericId(innerLoopNum) + elseif type == :chargePolar + diagid = PolarId(innerLoopNum, ChargeCharge, k=extK, t=tau_labels[extIndex]) + elseif type == :spinPolar + diagid = PolarId(innerLoopNum, SpinSpin, k=extK, t=tau_labels[extIndex]) + elseif type == :sigma + diagid = SigmaId(innerLoopNum, isDynamic ? Dynamic : Instant, k=extK, t=tau_labels[extIndex]) + else + end + if isempty(interactions) + # no interaction, so there should be only one spin factor. + g = Graph([fermi_greenProd], subgraph_factors=spinfactors_existed * symfactor, operator=IR.Sum(), properties=diagid) + else + inters = Graph(interactions, subgraph_factors=spinfactors_existed * symfactor, operator=IR.Sum()) + g = IR.multi_product(fermi_greenProd, inters, properties=diagid) + end + + return g +end + function read_onediagram!(io::IO, GNum::Int, verNum::Int, loopNum::Int, extIndex::Vector{Int}, labelProd::LabelProduct, spinPolarPara::Float64=0.0; diagType=:polar, maxLoopNum::Int=loopNum, splitter="|", offset::Int=-1, offset_ver4::Int=0, staticBose::Bool=true) From 3e85ec2560bed20a0a80634c8aaa42ab029a8e86 Mon Sep 17 00:00:00 2001 From: houpc Date: Fri, 19 Apr 2024 17:29:21 +0800 Subject: [PATCH 2/3] fix diagsGV() conflict --- src/frontend/GV.jl | 8 ++++---- 1 file changed, 4 insertions(+), 4 deletions(-) diff --git a/src/frontend/GV.jl b/src/frontend/GV.jl index 762409ef..2779f72e 100644 --- a/src/frontend/GV.jl +++ b/src/frontend/GV.jl @@ -37,8 +37,8 @@ include("GV_diagrams/readfile.jl") # Arguments: - `type` (Symbol): The type of the diagrams, including `:spinPolar`, `:chargePolar`, `:sigma`, `:green`, or `:freeEnergy`. - `order` (Int): The order of the diagrams without counterterms. -- `GOrder` (Int, optional): The order of self-energy counterterms (defaults to 0). -- `VerOrder` (Int, optional): The order of interaction counterterms (defaults to 0). +- `GOrder` (Int): The order of self-energy counterterms. +- `VerOrder` (Int): The order of interaction counterterms. - `labelProd` (Union{Nothing,LabelProduct}=nothing, optional): The initial cartesian QuantumOperator.label product (defaults to `nothing`). - `spinPolarPara` (Float64, optional): The spin-polarization parameter (n_up - n_down) / (n_up + n_down) (defaults to `0.0`). - `tau_labels`(Union{Nothing, Vector{Int}}, optional): The labels for the discrete time of each vertex. (defaults to `nothing`). @@ -49,9 +49,9 @@ A tuple `(diagrams, fermi_labelProd, bose_labelProd, extT_labels)` where - `labelProd` is a `LabelProduct` object containing the labels for the leaves of graphs, - `extT_labels` is a `Vector{Union{Tuple,Vector{Int}}}` object containing the external tau labels for each `FeynmanGraph` in `diagrams`. """ -function diagsGV(type::Symbol, order::Int, GOrder::Int=0, VerOrder::Int=0; +function diagsGV(type::Symbol, order::Int, GOrder::Int, VerOrder::Int; labelProd::Union{Nothing,LabelProduct}=nothing, spinPolarPara::Float64=0.0, - tau_labels::Union{Nothing,Vector{Int}}=nothing, filter::Vector{Filter}=[NoHartree]) + tau_labels::Union{Nothing,Vector{Int}}=nothing) if type == :spinPolar filename = string(@__DIR__, "/GV_diagrams/groups_spin/Polar$(order)_$(VerOrder)_$(GOrder).diag") elseif type == :chargePolar From 0455f1ac1d44c7587abbc0be07c16b9ce29e3c98 Mon Sep 17 00:00:00 2001 From: houpc Date: Fri, 19 Apr 2024 21:32:49 +0800 Subject: [PATCH 3/3] =?UTF-8?q?=E2=80=98?= MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit --- src/frontend/GV_diagrams/readfile.jl | 3 +++ 1 file changed, 3 insertions(+) diff --git a/src/frontend/GV_diagrams/readfile.jl b/src/frontend/GV_diagrams/readfile.jl index f7fbcef1..f31b30b6 100644 --- a/src/frontend/GV_diagrams/readfile.jl +++ b/src/frontend/GV_diagrams/readfile.jl @@ -521,6 +521,9 @@ function read_one_diagram!(type::Symbol, io::IO, GNum::Int, verNum::Int, loopNum spinFactors = _StringtoIntVector(readline(io)) extIndex = extIndex .- offset + if type == :sigma + extIndex[2] = findfirst(isequal(extIndex[1]), permutation) + end extNum = length(extIndex) extK = zeros(loopNum)