BP gauge for graph tensor network state#17
Conversation
|
Another problem: numerical instablity is observed, see #15 |
GiggleLiu
left a comment
There was a problem hiding this comment.
Very good. If you could add more comments or docstrings in the code, that would be perfect.
| @@ -0,0 +1,89 @@ | |||
| # the tensor network ansatz, |ket> | |||
| struct TensorNetworkAnsatz{TA, TB} | |||
| end | ||
| end | ||
|
|
||
| struct BPStep |
There was a problem hiding this comment.
Please add more descriptions to data structures to improved code readability.
| g::SimpleGraph{Int} # g is used to store the connections between the tensors | ||
| site_tensors::Vector{TA} # the tensors are stored as vector of arrays, corresponding to vertices of the graph. The virtual bounds are in order of its neighbors, the last dimension is the open dimension | ||
| gauge_tensors::Vector{TB} # corresponding to edges of the graph, listed in the order of edges(g), and enforce e.src < e.dst, real valued diagonal matrices | ||
| gauge_tensors_map::Dict{Tuple{Int, Int}, Int} # maps between the edge and the index of the gauge tensor |
There was a problem hiding this comment.
From tutorial perspective, I think maybe using SimpleEdge instead of Tuple{Int, Int} is easier for others to understand.
| return g | ||
| end | ||
|
|
||
| function chain(n::Int) |
|
|
||
| # codes about constructing the einsum expression | ||
|
|
||
| function all_eins(g::SimpleGraph{Int}, count::Int) |
There was a problem hiding this comment.
Bad design. Better not use count, which is a contextual information, as an argument. When I have a similar need of generating new indices, I usually pass an index store object to this function, e.g. https://github.com/CodingThrust/SimpleTDVP.jl/blob/fb946b589dc094b02587154d96efc8bc1ea99ec6/src/utils.jl#L1
| # 2:d + 1 are indices of the virtual indices of T | ||
| # d + 2:2d + 1 are indices of the virtual indices of T* | ||
| d = length(nebis) | ||
| T_ids = [2:d+1..., 1] |
There was a problem hiding this comment.
Consider using index store, which is concepturely easier.
| for i in eachindex(D) | ||
| (abs(D[i]) < 10 * eps(T)) && (D[i] = 10 * eps(T)) | ||
| end |
There was a problem hiding this comment.
| for i in eachindex(D) | |
| (abs(D[i]) < 10 * eps(T)) && (D[i] = 10 * eps(T)) | |
| end | |
| D .= max.(D, 10 * eps(T)) |
| # need to consider the pseudoinverse if D contains 0 | ||
| sqrt_D = sqrt.(D) | ||
| sqrt_M = U * diagm(sqrt_D) | ||
| sqrt_M_inv = diagm([!iszero(d) ? inv(d) : zero(typeof(d)) for d in sqrt_D]) * adjoint(U) |
There was a problem hiding this comment.
How can we have zero d? Also, zero(typeof(d)) is equivalent to zero(d).
There was a problem hiding this comment.
It is possible, since in some cases the ansatz is low ranked. For example, the boundary site of a chain with d_virtual > d_physical.
There was a problem hiding this comment.
But you have line 120-122. It will remove zero entries.
There was a problem hiding this comment.
Oh I forgot about that. 127 acutally is a pesudoinverse, but it does not work. I will remove that.
|
FYI: The checklist of open source software: openjournals/joss-reviews#5700 (comment) |
|
Revised according to most of the comments. |
What is done:
What is not done: