BMBPT Diagram¶
Routines and class for Bogoliubov MBPT diagrams.
- class adg.bmbpt.BmbptFeynmanDiagram(nx_graph, tag_num)[source]¶
Bases:
Diagram
Describes a BMBPT Feynman diagram with its related properties.
- Parameters
nx_graph (NetworkX MultiDiGraph) – The graph of interest.
tag_num (int) – The tag number associated to the diagram.
- attribute_expressions(time_diag)[source]¶
Attribute the correct Feynman and Goldstone expressions.
- Parameters
time_diag (TimeStructureDiagram) – The associated TSD.
- degrees¶
The ascendingly sorted degrees of the graph vertices.
- Type
tuple
- diag_exp¶
The Goldstone expression associated to the diagram.
- Type
str
- equivalent_permutations()[source]¶
Return the permutations generating equivalent diagrams.
- Returns
Vertices permutations as dictionnaries.
- Return type
list
- extract_integral()[source]¶
Return the integral part of the Feynman expression of the diag.
- Returns
The integral part of its Feynman expression.
- Return type
str
- extract_numerator()[source]¶
Return the numerator associated to a BMBPT graph.
- Returns
The numerator of the graph.
- Return type
str
- feynman_exp¶
The Feynman expression associated to the diagram.
- Type
str
- graph¶
The actual graph.
- Type
NetworkX MultiDiGraph
- has_crossing_sign()[source]¶
Return True for a minus sign associated with crossing propagators.
Use the fact that all lines propagate upwards and the canonical representation of the diagrams and vertices.
- Returns
- Encode for the sign factor associated with crossing
propagators.
- Return type
bool
- has_sign_factor()[source]¶
Return True if a sign factor is associated to the diagram.
Wrapper allowing for easy refactoring of expression code.
- Returns
The presence of a sign factor.
- Return type
bool
- hf_type¶
The Hartree-Fock, non-Hartree-Fock or Hartree-Fock for the energy operator only character of the graph.
- Type
str
- io_degrees¶
The sorted version of unsort_io_degrees.
- Type
tuple
- max_degree¶
The maximal degree of a vertex in the graph.
- Type
int
- multiplicity_symmetry_factor()[source]¶
Return the symmetry factor associated with propagators multiplicity.
- Returns
The symmetry factor associated with equivalent lines.
- Return type
str
- symmetry_factor()[source]¶
Return the overall symmetry factor of the diagram.
- Returns
The combination of all symmetry factors.
- Return type
str
- tags¶
The tag numbers associated to a diagram.
- Type
list
- time_tag¶
The tag number associated to the diagram’s associated TSD.
- Type
int
- time_tree_denominator(time_graph)[source]¶
Return the denominator for a time-tree graph.
- Parameters
time_graph (NetworkX MultiDiGraph) – Its associated time-structure graph.
- Returns
The denominator of the graph.
- Return type
str
- tsd_is_tree¶
The tree or non-tree character of the associated TSD.
- Type
bool
- two_or_three_body¶
The 2 or 3-body characted of the vertices.
- Type
int
- unique_id¶
A unique number associated to the diagram.
- Type
int
- unsort_degrees¶
The degrees of the graph vertices
- Type
tuple
- unsort_io_degrees¶
The list of in- and out-degrees for each vertex of the graph, stored in a (in, out) tuple.
- Type
tuple
- vert_exp¶
The expression associated to the vertices.
- Type
list
- property vertex_exchange_sym_factor¶
Return the symmetry factor associated with vertex exchange.
- Returns
The symmetry factor for vertex exchange.
- Return type
int
- vertex_expression(vertex)[source]¶
Return the expression associated to a given vertex.
- Parameters
vertex (int) – The vertex of interest in the graph.
- Returns
The LaTeX expression associated to the vertex.
- Return type
str
- write_diag_exps(latex_file, norder)[source]¶
Write the expressions associated to a diagram in the LaTeX file.
- Parameters
latex_file (file) – The LaTeX outputfile of the program.
norder (int) – The order in BMBPT formalism.
- write_graph(latex_file, directory, write_time)[source]¶
Write the BMBPT graph and its associated TSD to the LaTeX file.
- Parameters
latex_file (file) – The LaTeX output file of the program.
directory (str) – The path to the result folder.
write_time (bool) –
True
if we want informations on the associated TSDs.
- write_section(result, commands, section_flags)[source]¶
Write section and subsections for BMBPT result file.
- Parameters
result (file) – The LaTeX output file of the program.
commands (dict) – The flags associated with run management.
section_flags (dict) – UniqueIDs of diags starting each section.
- write_tsd_info(diagrams_time, latex_file)[source]¶
Write info related to the BMBPT associated TSD to the LaTeX file.
- Parameters
diagrams_time (list) – The associated TSDs.
latex_file (file) – The LaTeX output file of the program.
- write_vertices_values(latex_file, mapping)[source]¶
Write the qp energies associated to each vertex of the diag.
- Parameters
latex_file (file) – The LaTeX output file of the program.
mapping (dict) – A mapping between the vertices in the diagram and the vertices in its euivalent TSD, since permutations between vertices are possible.
- adg.bmbpt.check_topologically_equivalent(matrices, max_vertex)[source]¶
Exclude matrices that would spawn topologically equivalent graphs.
- Parameters
matrices (list) – Adjacency matrices to be checked.
max_vertex (int) – The maximum vertex which have been filled.
- Returns
The topologically unique matrices.
- Return type
list
>>> import numpy >>> mats = [numpy.array([[0, 2, 0, 0], [2, 0, 0, 1], [0, 0, 0, 0], [0, 0, 0, 0]]), numpy.array([[0, 0, 2, 0], [0, 0, 0, 0], [2, 0, 0, 1], [0, 0, 0, 0]])] >>> >>> mats = check_topologically_equivalent(mats, 2) >>> mats [array([[0, 2, 0, 0], [2, 0, 0, 1], [0, 0, 0, 0], [0, 0, 0, 0]])] >>> >>> mats = check_topologically_equivalent([], 2) >>> mats []
- adg.bmbpt.check_unconnected_spawn(matrices, max_filled_vertex)[source]¶
Exclude some matrices that would spawn unconnected diagrams.
Do several permutations among the rows and columns corresponding to already filled vertices, and check if one obtains a block-diagonal organisation, where the off-diagonals blocks connecting the already-filled and yet-unfilled parts of the matrix would be empty. In that case, remove the matrix.
- Parameters
matrices (list) – The adjacency matrices to be checked.
max_filled_vertex (int) – The furthest vertex until which the matrices have been filled.
>>> import numpy >>> mats = [numpy.array([[0, 2, 0], [2, 0, 0], [0, 0, 0]]), numpy.array([[0, 2, 1], [2, 0, 1], [0, 0, 0]])] >>> >>> check_unconnected_spawn(mats, 1) >>> mats [array([[0, 2, 1], [2, 0, 1], [0, 0, 0]])]
- adg.bmbpt.diagrams_generation(p_order, three_body_use, nbody_obs, canonical)[source]¶
Generate diagrams for BMBPT from bottom up.
- Parameters
p_order (int) – The BMBPT perturbative order of the studied diagrams.
three_body_use (bool) – Flag for the use of three-body forces.
nbody_obs (int) – N-body character of the obervable of interest.
canonical (bool) –
True
if one draws only canonical diagrams.
- Returns
NumPy arrays encoding the adjacency matrices of the graphs.
- Return type
list
>>> diags = diagrams_generation(1, False, 2, False) >>> len(diags) 2 >>> any(np.array_equal([[0, 4], [0, 0]], diag) for diag in diags) True >>> any(np.array_equal([[0, 2], [0, 0]], diag) for diag in diags) True >>> diags = diagrams_generation(1, True, 3, False) >>> len(diags) 3 >>> any(np.array_equal([[0, 6], [0, 0]], diag) for diag in diags) True >>> any(np.array_equal([[0, 4], [0, 0]], diag) for diag in diags) True >>> any(np.array_equal([[0, 2], [0, 0]], diag) for diag in diags) True >>> diags = diagrams_generation(2, False, 2, True) >>> len(diags) 2 >>> any(np.array_equal([[0, 2, 2], [0, 0, 2], [0, 0, 0]], d) for d in diags) True >>> any(np.array_equal([[0, 1, 1], [0, 0, 3], [0, 0, 0]], d) for d in diags) True
- adg.bmbpt.order_and_remove_topologically_equiv(matrices, max_vertex)[source]¶
Order the matrices in sub-list and remove topologically equivalent ones.
- Parameters
matrices (list) – The adjacency matrices to be checked.
max_vertex (int) – The maximum vertex which has been filled.
- Returns
The ordered topologically unique matrices.
- Return type
list
- adg.bmbpt.order_diagrams(diagrams)[source]¶
Order the BMBPT diagrams and return number of diags for each type.
- Parameters
diagrams (list) – Possibly redundant BmbptFeynmanDiagrams.
- Returns
- First element is the list of topologically unique, ordered
diagrams. Second element is a dict with the number of diagrams for each major type. Third element is a dict with the identifiers of diagrams starting each output file section.
- Return type
tuple
- adg.bmbpt.produce_expressions(diagrams, diagrams_time)[source]¶
Produce and store the expressions associated to the BMBPT diagrams.
- Parameters
diagrams (list) – The list of all BmbptFeynmanDiagrams.
diagrams_time (list) – Their associates TSDs.