关于磁矩特征的写入

以 Ti2O3 为示例,cif 文件为:

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
# generated using pymatgen
data_Ti2O3
_symmetry_space_group_name_H-M 'P 1'
_cell_length_a 5.16780216
_cell_length_b 5.16780216
_cell_length_c 9.60760300
_cell_angle_alpha 90.00000000
_cell_angle_beta 90.00000000
_cell_angle_gamma 120.00000000
_symmetry_Int_Tables_number 1
_chemical_formula_structural Ti2O3
_chemical_formula_sum 'Ti8 O12'
_cell_volume 222.20684815
_cell_formula_units_Z 4
loop_
_symmetry_equiv_pos_site_id
_symmetry_equiv_pos_as_xyz
1 'x, y, z'
loop_
_atom_site_type_symbol
_atom_site_label
_atom_site_symmetry_multiplicity
_atom_site_fract_x
_atom_site_fract_y
_atom_site_fract_z
_atom_site_occupancy
Ti Ti0 1 0.33333333 0.66666667 0.97654800 1
Ti Ti1 1 0.33333333 0.66666667 0.25000000 1
Ti Ti2 1 0.00000000 0.00000000 0.00000000 1
Ti Ti3 1 0.00000000 0.00000000 0.50000000 1
Ti Ti4 1 0.33333333 0.66666667 0.52345200 1
Ti Ti5 1 0.66666667 0.33333333 0.75000000 1
Ti Ti6 1 0.66666667 0.33333333 0.02345200 1
Ti Ti7 1 0.66666667 0.33333333 0.47654800 1
O O8 1 0.02585733 0.67555467 0.11641700 1
O O9 1 0.97414267 0.64969733 0.61641700 1
O O10 1 0.02585733 0.35030267 0.38358300 1
O O11 1 0.97414267 0.32444533 0.88358300 1
O O12 1 0.32444533 0.35030267 0.11641700 1
O O13 1 0.35030267 0.32444533 0.61641700 1
O O14 1 0.32444533 0.97414267 0.38358300 1
O O15 1 0.35030267 0.02585733 0.88358300 1
O O16 1 0.64969733 0.67555467 0.38358300 1
O O17 1 0.67555467 0.64969733 0.88358300 1
O O18 1 0.64969733 0.97414267 0.11641700 1
O O19 1 0.67555467 0.02585733 0.61641700 1

可以在下面的原子 loop 矩阵中加入每个原子对应的标量磁矩 _atom_site_moment

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
# 头部信息不变
loop_
_atom_site_type_symbol
_atom_site_label
_atom_site_symmetry_multiplicity
_atom_site_fract_x
_atom_site_fract_y
_atom_site_fract_z
_atom_site_occupancy
_atom_site_moment
Ti Ti0 1 0.33333333 0.66666667 0.97654800 1 5
Ti Ti1 1 0.33333333 0.66666667 0.25000000 1 5
Ti Ti2 1 0.00000000 0.00000000 0.00000000 1 2
Ti Ti3 1 0.00000000 0.00000000 0.50000000 1 2
Ti Ti4 1 0.33333333 0.66666667 0.52345200 1 2
Ti Ti5 1 0.66666667 0.33333333 0.75000000 1 2
Ti Ti6 1 0.66666667 0.33333333 0.02345200 1 1
Ti Ti7 1 0.66666667 0.33333333 0.47654800 1 1
O O8 1 0.02585733 0.67555467 0.11641700 1 0
O O9 1 0.97414267 0.64969733 0.61641700 1 0
O O10 1 0.02585733 0.35030267 0.38358300 1 0
O O11 1 0.97414267 0.32444533 0.88358300 1 0
O O12 1 0.32444533 0.35030267 0.11641700 1 0
O O13 1 0.35030267 0.32444533 0.61641700 1 0
O O14 1 0.32444533 0.97414267 0.38358300 1 0
O O15 1 0.35030267 0.02585733 0.88358300 1 0
O O16 1 0.64969733 0.67555467 0.38358300 1 0
O O17 1 0.67555467 0.64969733 0.88358300 1 0
O O18 1 0.64969733 0.97414267 0.11641700 1 0
O O19 1 0.67555467 0.02585733 0.61641700 1 0

或者还可以写入矢量磁矩,但是需要 3 个分量

1
2
3
4
5
loop_
_atom_site_label
_atom_site_moment_Cartn_x
_atom_site_moment_Cartn_y
_atom_site_moment_Cartn_z

如何获取每个原子的磁矩信息?

我查阅https://next-gen.materialsproject.org/materials/mp-776655?formula=Ti2O3#more

中含有磁矩信息

获取 Ti2O3 全部信息的方式为

1
2
3
from mp_api.client import MPRester
with MPRester(api_key="igRHy7zYOKzWD18jY76XtjbwEpRl6SoH") as mpr:
data = mpr.materials.search(material_ids=["mp-776655"])

获取之后打印全部 data 信息得到

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
[MPDataDoc<MaterialsDoc>(
builder_meta=EmmetMeta(emmet_version='0.84.3rc4',
pymatgen_version='2024.11.13',
run_id='ccba9c49-e919-43c2-8d4f-35b2d60eb94b',
batch_id=None, database_version='2025.06.09',
build_date=datetime.datetime(2024, 11, 21, 22, 57, 42, 12000,
tzinfo=datetime.timezone.utc),
license='BY-C'),
nsites=20,
elements=[Element O, Element Ti],
nelements=2,
composition=Composition('Ti8 O12'),
composition_reduced=Composition('Ti2 O3'),
formula_pretty='Ti2O3',
formula_anonymous='A2B3',
chemsys='O-Ti',
volume=217.75988736026858,
density=4.384148726776235,
density_atomic=10.887994368013429,
symmetry=SymmetryData(crystal_system=<CrystalSystem.trig: 'Trigonal'>, symbol='P-31c', number=163, point_group='-3m', symprec=0.1, angle_tolerance=5.0, version='2.5.0'),
material_id=MPID(mp-776655),
structure=Structure Summary
Lattice
abc : 5.151006976476778 5.150980372462452 9.476884040021895
angles : 90.0000031964729 90.0002147954698 120.00007956211232
volume : 217.75988736026858
A : 2.57550686 -4.46090095 -9.72e-06
B : 2.57548062 4.46088538 -2.9e-07
C : -1.764e-05 1.019e-05 9.47688404
pbc : True True True
PeriodicSite: Ti (2.575, 1.487, 9.24) [0.3333, 0.6667, 0.975]
PeriodicSite: Ti (2.575, 1.487, 2.369) [0.3333, 0.6667, 0.25]
PeriodicSite: Ti (-1.857e-05, 1.073e-05, 9.477) [-2.4e-07, -1.2e-07, 1.0]
PeriodicSite: Ti (-4.39e-06, 1.518e-05, 4.739) [-2.7e-07, 1.99e-06, 0.5]
PeriodicSite: Ti (2.575, 1.487, 4.975) [0.3333, 0.6667, 0.525]
PeriodicSite: Ti (2.575, -1.487, 7.108) [0.6667, 0.3333, 0.75]
PeriodicSite: Ti (2.576, -1.487, 0.2363) [0.6667, 0.3333, 0.02494]
PeriodicSite: Ti (2.576, -1.487, 4.502) [0.6667, 0.3333, 0.4751]
PeriodicSite: O (1.803, 2.894, 1.104) [0.02568, 0.6744, 0.1165]
PeriodicSite: O (4.18, -1.453, 5.843) [0.9743, 0.6487, 0.6166]
PeriodicSite: O (0.971, 1.452, 3.634) [0.0257, 0.3513, 0.3835]
PeriodicSite: O (3.348, -2.894, 8.372) [0.9743, 0.3256, 0.8834]
PeriodicSite: O (1.743, 0.1145, 1.104) [0.3256, 0.3513, 0.1165]
PeriodicSite: O (1.743, -0.1146, 5.843) [0.3513, 0.3256, 0.6166]
PeriodicSite: O (3.348, 2.894, 3.634) [0.3256, 0.9743, 0.3835]
PeriodicSite: O (0.9709, -1.452, 8.372) [0.3513, 0.02569, 0.8834]
PeriodicSite: O (3.408, 0.1146, 3.634) [0.6487, 0.6744, 0.3835]
PeriodicSite: O (3.408, -0.1146, 8.372) [0.6744, 0.6487, 0.8834]
PeriodicSite: O (4.18, 1.452, 1.104) [0.6487, 0.9743, 0.1165]
PeriodicSite: O (1.803, -2.894, 5.843) [0.6744, 0.02568, 0.6166],
deprecated=False,
deprecation_reasons=None,
initial_structures=[Structure Summary
Lattice
abc : 5.167048146727587 5.167048146727587 9.607603
angles : 90.0 90.0 120.00000187871662
volume : 222.14200594619678
A : 2.583523999999999 -4.474795 0.0
B : 2.583523999999999 4.474795 0.0
C : 0.0 0.0 9.607603
pbc : True True True
PeriodicSite: Ti (2.584, 1.492, 9.382) [0.3333, 0.6667, 0.9765]
PeriodicSite: Ti (2.584, 1.492, 2.402) [0.3333, 0.6667, 0.25]
PeriodicSite: Ti (0.0, 0.0, 0.0) [0.0, 0.0, 0.0]
PeriodicSite: Ti (0.0, 0.0, 4.804) [0.0, 0.0, 0.5]
PeriodicSite: Ti (2.584, 1.492, 5.029) [0.3333, 0.6667, 0.5235]
PeriodicSite: Ti (2.584, -1.492, 7.206) [0.6667, 0.3333, 0.75]
PeriodicSite: Ti (2.584, -1.492, 0.2253) [0.6667, 0.3333, 0.02345]
PeriodicSite: Ti (2.584, -1.492, 4.578) [0.6667, 0.3333, 0.4765]
PeriodicSite: O (1.743, 0.1157, 1.118) [0.3244, 0.3503, 0.1164]
PeriodicSite: O (1.743, -0.1157, 5.922) [0.3503, 0.3244, 0.6164]
PeriodicSite: O (3.424, 0.1157, 3.685) [0.6497, 0.6756, 0.3836]
PeriodicSite: O (3.424, -0.1157, 8.489) [0.6756, 0.6497, 0.8836]
PeriodicSite: O (4.195, 1.452, 1.118) [0.6497, 0.9741, 0.1164]
PeriodicSite: O (1.812, -2.907, 5.922) [0.6756, 0.02586, 0.6164]
PeriodicSite: O (0.9718, 1.452, 3.685) [0.02586, 0.3503, 0.3836]
PeriodicSite: O (3.355, -2.907, 8.489) [0.9741, 0.3244, 0.8836]
PeriodicSite: O (3.355, 2.907, 3.685) [0.3244, 0.9741, 0.3836]
PeriodicSite: O (0.9718, -1.452, 8.489) [0.3503, 0.02586, 0.8836]
PeriodicSite: O (1.812, 2.907, 1.118) [0.02586, 0.6756, 0.1164]
PeriodicSite: O (4.195, -1.452, 5.922) [0.9741, 0.6497, 0.6164]],
task_ids=[MPID(mp-2739261), MPID(mp-2739212), MPID(mp-888770), MPID(mp-1784547), MPID(mp-902522), MPID(mp-903281), MPID(mp-1367622), MPID(mp-776655)],
deprecated_tasks=[],
calc_types={'mp-902522': <CalcType.GGA_NSCF_Line: 'GGA NSCF Line'>, 'mp-776655': <CalcType.GGA_Structure_Optimization: 'GGA Structure Optimization'>, 'mp-888770': <CalcType.GGA_Static: 'GGA Static'>, 'mp-1367622': <CalcType.GGA_Static: 'GGA Static'>, 'mp-1784547': <CalcType.GGA_Static: 'GGA Static'>, 'mp-903281': <CalcType.GGA_NSCF_Line: 'GGA NSCF Line'>, 'mp-2739212': <CalcType.PBEsol_Structure_Optimization: 'PBEsol Structure Optimization'>, 'mp-2739261': <CalcType.r2SCAN_Structure_Optimization: 'r2SCAN Structure Optimization'>},
last_updated=datetime.datetime(2021, 7, 25, 14, 16, 23, 198000, tzinfo=datetime.timezone.utc),
created_at=datetime.datetime(2014, 2, 22, 5, 1, 31, tzinfo=datetime.timezone.utc),
origins=[PropertyOrigin(name='structure', task_id=MPID(mp-2739261), last_updated=datetime.datetime(2021, 7, 25, 14, 16, 23, 198000, tzinfo=datetime.timezone.utc))],
warnings=[],
task_types={'mp-902522': <TaskType.NSCF_Line: 'NSCF Line'>, 'mp-776655': <TaskType.Structure_Optimization: 'Structure Optimization'>, 'mp-888770': <TaskType.Static: 'Static'>, 'mp-1367622': <TaskType.Static: 'Static'>, 'mp-1784547': <TaskType.Static: 'Static'>, 'mp-903281': <TaskType.NSCF_Line: 'NSCF Line'>, 'mp-2739212': <TaskType.Structure_Optimization: 'Structure Optimization'>, 'mp-2739261': <TaskType.Structure_Optimization: 'Structure Optimization'>},
run_types={'mp-902522': <RunType.GGA: 'GGA'>, 'mp-776655': <RunType.GGA: 'GGA'>, 'mp-888770': <RunType.GGA: 'GGA'>, 'mp-1367622': <RunType.GGA: 'GGA'>, 'mp-1784547': <RunType.GGA: 'GGA'>, 'mp-903281': <RunType.GGA: 'GGA'>, 'mp-2739212': <RunType.PBEsol: 'PBEsol'>, 'mp-2739261': <RunType.r2SCAN: 'r2SCAN'>},
entries=BlessedCalcs(GGA=mp-776655-GGA ComputedStructureEntry - Ti8 O12 (Ti2O3)
Energy (Uncorrected) = -178.7848 eV (-8.9392 eV/atom)
Correction = 0.0000 eV (0.0000 eV/atom)
Energy (Final) = -178.7848 eV (-8.9392 eV/atom)
Energy Adjustments:
None
Parameters:
potcar_spec = [{'titel': 'PAW_PBE Ti_pv 07Sep2000', 'hash': '70bc3ea8bf68f10e7e1e4721bb91972a', 'summary_stats': None}, {'titel': 'PAW_PBE O 08Apr2002', 'hash': '7a25bc5b9a5393f46600a4939d357982', 'summary_stats': None}]
run_type = GGA
is_hubbard = False
hubbards = None
Data:
oxide_type = oxide
aspherical = True
last_updated = 2024-11-21 22:57:41.717606+00:00
task_id = mp-1784547
material_id = mp-776655, GGA_U=None, PBESol=mp-776655-PBEsol ComputedStructureEntry - Ti8 O12 (Ti2O3)
Energy (Uncorrected) = -186.9983 eV (-9.3499 eV/atom)
Correction = 0.0000 eV (0.0000 eV/atom)
Energy (Final) = -186.9983 eV (-9.3499 eV/atom)
Energy Adjustments:
None
Parameters:
potcar_spec = [{'titel': 'PAW_PBE Ti_pv 07Sep2000', 'hash': 'cdd047e254a1247dbf5c1bb0563a75b9', 'summary_stats': None}, {'titel': 'PAW_PBE O 08Apr2002', 'hash': '9bb4b91e6c47f70fd2bce603bd5d6832', 'summary_stats': None}]
run_type = PBEsol
is_hubbard = False
hubbards = None
Data:
oxide_type = oxide
aspherical = True
last_updated = 2024-11-21 22:57:41.767313+00:00
task_id = mp-2739212
material_id = mp-776655, SCAN=None, R2SCAN=mp-776655-r2SCAN ComputedStructureEntry - Ti8 O12 (Ti2O3)
Energy (Uncorrected) = -236.5932 eV (-11.8297 eV/atom)
Correction = 0.0000 eV (0.0000 eV/atom)
Energy (Final) = -236.5932 eV (-11.8297 eV/atom)
Energy Adjustments:
None
Parameters:
potcar_spec = [{'titel': 'PAW_PBE Ti_pv 07Sep2000', 'hash': 'cdd047e254a1247dbf5c1bb0563a75b9', 'summary_stats': None}, {'titel': 'PAW_PBE O 08Apr2002', 'hash': '9bb4b91e6c47f70fd2bce603bd5d6832', 'summary_stats': None}]
run_type = r2SCAN
is_hubbard = False
hubbards = None
Data:
oxide_type = oxide
aspherical = True
last_updated = 2024-11-21 22:57:41.802864+00:00
task_id = mp-2739261
material_id = mp-776655, HSE=None),
fields_not_requested=[]
)]

发现 structure 字段与我本地的 cif 文件不太一样,他多了一些额外的字段

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
loop_
_atom_site_type_symbol
_atom_site_label
_atom_site_symmetry_multiplicity
_atom_site_fract_x
_atom_site_fract_y
_atom_site_fract_z
_atom_site_occupancy
Ti Ti0 1 0.33333333 0.66666667 0.97654800 1
Ti Ti1 1 0.33333333 0.66666667 0.25000000 1
Ti Ti2 1 0.00000000 0.00000000 0.00000000 1
Ti Ti3 1 0.00000000 0.00000000 0.50000000 1
Ti Ti4 1 0.33333333 0.66666667 0.52345200 1
Ti Ti5 1 0.66666667 0.33333333 0.75000000 1
Ti Ti6 1 0.66666667 0.33333333 0.02345200 1
Ti Ti7 1 0.66666667 0.33333333 0.47654800 1
O O8 1 0.02585733 0.67555467 0.11641700 1
O O9 1 0.97414267 0.64969733 0.61641700 1
O O10 1 0.02585733 0.35030267 0.38358300 1
O O11 1 0.97414267 0.32444533 0.88358300 1
O O12 1 0.32444533 0.35030267 0.11641700 1
O O13 1 0.35030267 0.32444533 0.61641700 1
O O14 1 0.32444533 0.97414267 0.38358300 1
O O15 1 0.35030267 0.02585733 0.88358300 1
O O16 1 0.64969733 0.67555467 0.38358300 1
O O17 1 0.67555467 0.64969733 0.88358300 1
O O18 1 0.64969733 0.97414267 0.11641700 1
O O19 1 0.67555467 0.02585733 0.61641700 1

于是我就打印下来看看

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
Full Formula (Ti8 O12)
Reduced Formula: Ti2O3
abc : 5.151007 5.150980 9.476884
angles: 90.000003 90.000215 120.000080
pbc : True True True
Sites (20)
# SP a b c magmom
--- ---- --------- --------- -------- --------
0 Ti 0.333332 0.666663 0.975046 0.514
1 Ti 0.333334 0.666667 0.249997 -0.112
2 Ti -0 -0 0.999985 0.906
3 Ti -0 2e-06 0.500008 0.906
4 Ti 0.333322 0.666661 0.524963 0.514
5 Ti 0.666666 0.333333 0.750007 -0.112
6 Ti 0.666668 0.333337 0.024937 0.514
7 Ti 0.666677 0.333337 0.475051 0.514
8 O 0.025681 0.674403 0.116537 -0.007
9 O 0.974318 0.648705 0.616563 -0.007
10 O 0.0257 0.351301 0.383472 -0.007
11 O 0.974308 0.325583 0.883429 -0.007
12 O 0.325599 0.351278 0.116537 -0.007
13 O 0.351293 0.32561 0.61656 -0.007
14 O 0.3256 0.974301 0.383477 -0.007
15 O 0.351274 0.02569 0.88343 -0.007
16 O 0.648702 0.674402 0.383477 -0.007
17 O 0.674414 0.648724 0.883429 -0.007
18 O 0.648724 0.974321 0.116535 -0.007
19 O 0.674388 0.025681 0.61656 -0.007

可见 structure 中是含有原子磁矩的,于是我就将其重新写 cif,在这里我用了 2 种方法

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
# 方法一
structure.to("Ti2O3.cif")
# 方法二
from pymatgen.io.cif import CifWriter
cif_writer = CifWriter(structure, write_magmoms=True)

### Cifwriter的参数
class CifWriter:
"""A wrapper around CifFile to write CIF files from pymatgen Structure."""

def __init__(
self,
struct: Structure | IStructure,
symprec: float | None = None,
write_magmoms: bool = False,
significant_figures: int = 8,
angle_tolerance: float = 5,
refine_struct: bool = True,
write_site_properties: bool = False,
) -> None:
"""
Args:
struct (Structure): structure to write.
symprec (float): If not none, finds the symmetry of the structure
and writes the CIF with symmetry information. Passes symprec
to the SpacegroupAnalyzer. See also refine_struct.
write_magmoms (bool): If True, will write magCIF file. Incompatible
with symprec
significant_figures (int): Specifies precision for formatting of floats.
Defaults to 8.
angle_tolerance (float): Angle tolerance for symmetry finding. Passes
angle_tolerance to the SpacegroupAnalyzer. Used only if symprec
is not None.
refine_struct: Used only if symprec is not None. If True, get_refined_structure
is invoked to convert input structure from primitive to conventional.
write_site_properties (bool): Whether to write the Structure.site_properties
to the CIF as _atom_site_{property name}. Defaults to False.
=====================================================================
struct (Structure): 待写入的晶体结构。
symprec (float): 若不为None,则分析结构的对称性并写入带对称性信息的CIF文件。
该参数会传递给空间群分析器(SpacegroupAnalyzer)。另见refine_struct参数说明。
write_magmoms (bool): 若为True,则写入磁矩信息生成magCIF文件。与symprec参数不兼容。
significant_figles (int): 浮点数格式化精度,默认保留8位有效数字。
angle_tolerance (float): 对称性分析的角度容差参数,传递给SpacegroupAnalyzer。
仅在symprec不为None时生效。
refine_struct (bool): 仅在symprec不为None时生效。若为True,则调用
get_refined_structure方法将原结构从原胞转换为晶胞常规形式。
write_site_properties (bool): 是否将Structure.site_properties中的属性
写入CIF的_atom_site_{属性名}字段,默认为False。
"""

第一种方法是没有用的,第二种方法可以写入,得到以下矢量磁矩

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
loop_
_atom_site_moment_label
_atom_site_moment_crystalaxis_x
_atom_site_moment_crystalaxis_y
_atom_site_moment_crystalaxis_z
Ti0 0.00000128 0.00000064 0.51400000
Ti1 -0.00000028 -0.00000014 -0.11200000
Ti2 0.00000225 0.00000112 0.90600000
Ti3 0.00000225 0.00000112 0.90600000
Ti4 0.00000128 0.00000064 0.51400000
Ti5 -0.00000028 -0.00000014 -0.11200000
Ti6 0.00000128 0.00000064 0.51400000
Ti7 0.00000128 0.00000064 0.51400000
O8 -0.00000002 0.00000000 -0.00700000
O9 -0.00000002 0.00000000 -0.00700000
O10 -0.00000002 0.00000000 -0.00700000
O11 -0.00000002 0.00000000 -0.00700000
O12 -0.00000002 0.00000000 -0.00700000
O13 -0.00000002 0.00000000 -0.00700000
O14 -0.00000002 0.00000000 -0.00700000
O15 -0.00000002 0.00000000 -0.00700000
O16 -0.00000002 0.00000000 -0.00700000
O17 -0.00000002 0.00000000 -0.00700000
O18 -0.00000002 0.00000000 -0.00700000
O19 -0.00000002 0.00000000 -0.00700000

接下来尝试获取标量磁矩,试了很多方法,都失败了,而且很复杂

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
# 提取矢量磁矩数据
import numpy as np
vector_magmoms = structure.site_properties.get("magmom", [])

# 计算标量磁矩(模长)
scalar_magmoms = []
for magmom in vector_magmoms:
# 检查磁矩是否是矢量(列表/数组)还是标量
if isinstance(magmom, (list, tuple, np.ndarray)):
# 计算矢量模长作为标量磁矩
scalar_magmoms.append(np.linalg.norm(magmom))
else:
# 如果已经是标量,直接使用
scalar_magmoms.append(magmom)

# 使用CifWriter但添加自定义标量磁矩列
cif_writer = CifWriter(structure)
cif_file = cif_writer.cif_file
# 获取第一个数据块
data_block = next(iter(cif_file.data.values()))
# 添加标量磁矩数据
atom_labels = [f"Ti{i}" if s.species_string == "Ti" else f"O{i}"
for i, s in enumerate(structure)]

# 添加标量磁矩循环
data_block.add_loop([
["_atom_site_label"] + atom_labels,
["_atom_site_scalar_magmom"] + [f"{m:.6f}" for m in scalar_magmoms]
])

# 保存修改后的CIF
alt_filename = f"{MP_ID}_alt_scalar_magmoms.cif"
cif_file.write_file(alt_filename)
print(f"备选方法生成CIF文件: {alt_filename}")

# 打印验证信息
print("\n标量磁矩值:")
for i, magmom in enumerate(scalar_magmoms):
element = structure[i].species_string
print(f"原子 {i} ({element}): {magmom:.6f} μB")

如何读取这些磁矩并加入训练?

首先,我按照原来的方式读取了原始数据和加入了磁矩的数据

1
2
3
4
5
6
7
8
9
10
11
from csat.crystal_data import CIFData, crystal_graph_list
from csat.subgraph_data import KhopGraphDataset
data = CIFData("Ti2O3", target_name='is_Magnetic' )
print('这是读取cif文件后的数据')
graph = crystal_graph_list(data)
print('这是图数据')
print(graph[0])
print(graph[0].x)
graph_sub = KhopGraphDataset(graph, k_hop=2)
print('这是子图数据')
print(graph_sub[0])

两次结果如下

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
这是读取cif文件后的数据
100%|██████████| 1/1 [00:00<00:00, 29.05it/s]
这是图数据
Data(x=[20, 92], edge_index=[2, 196], edge_attr=[196, 41], y=[1], id='Ti2O3')
tensor([[0., 0., 0., ..., 0., 0., 0.],
[0., 0., 0., ..., 0., 0., 0.],
[0., 0., 0., ..., 0., 0., 0.],
...,
[0., 0., 0., ..., 0., 0., 0.],
[0., 0., 0., ..., 0., 0., 0.],
[0., 0., 0., ..., 0., 0., 0.]])
Extracting 2-hop subgraphs...
100%|██████████| 1/1 [00:00<00:00, 77.63it/s]
Done!
这是子图数据
Data(x=[20, 92], edge_index=[2, 196], edge_attr=[196, 41], y=[1, 1], id='Ti2O3', degree=[20], complete_edge_index=[2, 400], subgraph_edge_index=[2, 3756], num_subgraph_nodes=390, subgraph_node_idx=[390], subgraph_indicator=[390]
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
这是读取cif文件后的数据
100%|██████████| 1/1 [00:00<00:00, 27.14it/s]
这是图数据
Data(x=[20, 92], edge_index=[2, 196], edge_attr=[196, 41], y=[1], id='Ti2O3')
tensor([[0., 0., 0., ..., 0., 0., 0.],
[0., 0., 0., ..., 0., 0., 0.],
[0., 0., 0., ..., 0., 0., 0.],
...,
[0., 0., 0., ..., 0., 0., 0.],
[0., 0., 0., ..., 0., 0., 0.],
[0., 0., 0., ..., 0., 0., 0.]])
Extracting 2-hop subgraphs...
100%|██████████| 1/1 [00:00<00:00, 77.26it/s]
Done!
这是子图数据
Data(x=[20, 92], edge_index=[2, 196], edge_attr=[196, 41], y=[1, 1], id='Ti2O3', degree=[20], complete_edge_index=[2, 400], subgraph_edge_index=[2, 3756], num_subgraph_nodes=390, subgraph_node_idx=[390], subgraph_indicator=[390])

没有任何区别,说明没有读取到磁矩信息,于是我就要从读取方式上着手,追本溯源。

在原来项目文件中,读取 cif 结构的部分如下:

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
from pymatgen.core.structure import Structure

crystal = Structure.from_file(os.path.join(self.root_dir,
cif_id + '.cif'))
atom_fea = np.vstack([self.ari.get_atom_fea(crystal[i].specie.number)
for i in range(len(crystal))])
atom_fea = torch.Tensor(atom_fea)
all_nbrs = crystal.get_all_neighbors(self.radius, include_index=True)
all_nbrs = [sorted(nbrs, key=lambda x: x[1]) for nbrs in all_nbrs]
nbr_fea_idx, nbr_fea = [], []
for nbr in all_nbrs:
if len(nbr) < self.max_num_nbr:
warnings.warn('{} not find enough neighbors to build graph. '
'If it happens frequently, consider increase '
'radius.'.format(cif_id))
nbr_fea_idx.append(list(map(lambda x: x[2], nbr)) +
[0] * (self.max_num_nbr - len(nbr)))
nbr_fea.append(list(map(lambda x: x[1], nbr)) +
[self.radius + 1.] * (self.max_num_nbr -
len(nbr)))
else:
nbr_fea_idx.append(list(map(lambda x: x[2],
nbr[:self.max_num_nbr])))
nbr_fea.append(list(map(lambda x: x[1],
nbr[:self.max_num_nbr])))
nbr_fea_idx, nbr_fea = np.array(nbr_fea_idx), np.array(nbr_fea)
nbr_fea = self.gdf.expand(nbr_fea)
atom_fea = torch.Tensor(atom_fea)
nbr_fea = torch.Tensor(nbr_fea)
nbr_fea_idx = torch.LongTensor(nbr_fea_idx)
target = torch.Tensor([float(target_value)])
return (atom_fea, nbr_fea, nbr_fea_idx), target, cif_id

我们来一步一步测试

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
from pymatgen.core.structure import Structure
crystal = Structure.from_file('Ti2O3_m/Ti2O3.cif')
print(crystal)

'''
Full Formula (Ti8 O12)
Reduced Formula: Ti2O3
abc : 5.151007 5.150980 9.476884
angles: 90.000003 90.000215 120.000080
pbc : True True True
Sites (20)
# SP a b c magmom
--- ---- -------- -------- -------- --------
0 Ti 0.333333 0.666667 0.975046 0.514
1 Ti 0.333333 0.666667 0.249997 -0.112
2 Ti 1 1 0.999985 0.906
3 Ti 1 2e-06 0.500008 0.906
4 Ti 0.333333 0.666667 0.524963 0.514
5 Ti 0.666667 0.333333 0.750007 -0.112
6 Ti 0.666667 0.333333 0.024937 0.514
7 Ti 0.666667 0.333333 0.475051 0.514
8 O 0.025681 0.674403 0.116537 -0.007
9 O 0.974318 0.648705 0.616563 -0.007
10 O 0.0257 0.351301 0.383472 -0.007
11 O 0.974308 0.325583 0.883429 -0.007
12 O 0.325599 0.351278 0.116537 -0.007
13 O 0.351293 0.32561 0.61656 -0.007
14 O 0.3256 0.974301 0.383477 -0.007
15 O 0.351274 0.02569 0.88343 -0.007
16 O 0.648702 0.674402 0.383477 -0.007
17 O 0.674414 0.648724 0.883429 -0.007
18 O 0.648724 0.974321 0.116535 -0.007
19 O 0.674388 0.025681 0.61656 -0.007
'''
这里是利用 pymatgen 读取 structure,可以看见跟之前获取时候打印的是一模一样的,而且是含有磁矩的。
接下来按照谢天的方法来读取原子特征
```python
import numpy as np
import json
import torch
# 原子特征初始化器
# 原子特征初始化器
class AtomInitializer(object):
"""
Base class for intializing the vector representation for atoms.
该类的作用是初始化原子的特征表示。
初始化过程是通过一个字典将每种元素的原子特征表示存储在内存中。
"""

def __init__(self, atom_types):
self.atom_types = set(atom_types) # 存储原子类型的集合
self._embedding = {} # # 核心存储:原子类型 -> 特征向量

def get_atom_fea(self, atom_type):
assert atom_type in self.atom_types
return self._embedding[atom_type] # 原子序数作为键获取特征向量

def load_state_dict(self, state_dict):
self._embedding = state_dict
self.atom_types = set(self._embedding.keys())
self._decodedict = {idx: atom_type for atom_type, idx in
self._embedding.items()}

def state_dict(self):
return self._embedding

def decode(self, idx):
if not hasattr(self, '_decodedict'):
self._decodedict = {idx: atom_type for atom_type, idx in
self._embedding.items()}
return self._decodedict[idx]

class AtomCustomJSONInitializer(AtomInitializer):
"""
继承自 AtomInitializer,专门通过一个 JSON 文件加载原子特征,该文件包含了每种元素的特征向量。
"""

def __init__(self, elem_embedding_file):
with open(elem_embedding_file) as f:
elem_embedding = json.load(f) # 读取 JSON 文件

# 将 JSON 中的键转换为整数类型,并将值转换为 NumPy 数组
# 这样可以确保原子序数作为键,特征向量作为值
elem_embedding = {int(key): value for key, value
in elem_embedding.items()}

# 初始化父类,传入原子类型集合
# 这里的原子类型是 JSON 文件中的键(整数形式的原子序数)
# 例如:{6: [0.1, 0.2, ... ], 8: [0.3, 0.4, ...]}
atom_types = set(elem_embedding.keys())
# 调用父类构造函数
super(AtomCustomJSONInitializer, self).__init__(atom_types)
# 将 JSON 中的特征向量存储到 _embedding 字典中
# 这里的键是整数形式的原子序数,值是对应的特征向量
# 例如:{6: np.array([0.1, 0.2
for key, value in elem_embedding.items():
self._embedding[key] = np.array(value, dtype=float)


ari = AtomCustomJSONInitializer('Ti2O3_m/atom_init.json')
# 遍历晶体中的每个原子
atom_fea = np.vstack([ari.get_atom_fea(crystal[i].specie.number) # 获取原子序数的特征向量
for i in range(len(crystal))])
atom_fea = torch.Tensor(atom_fea)
print('这是原子特征')
print(atom_fea)

'''结果
这是原子特征
tensor([[0., 0., 0., ..., 0., 0., 0.],
[0., 0., 0., ..., 0., 0., 0.],
[0., 0., 0., ..., 0., 0., 0.],
...,
[0., 0., 0., ..., 0., 0., 0.],
[0., 0., 0., ..., 0., 0., 0.],
[0., 0., 0., ..., 0., 0., 0.]])
'''

从这里就可以发现,原子特征是按照 atom_init.json 进行编码的。

这段代码的核心功能是:将晶体结构中的每个原子转化为数值化特征向量,最终生成一个PyTorch张量。

具体流程:

  • 从JSON文件加载预定义的原子特征向量(每个原子序数对应一个特征向量)
  • 遍历晶体中的每个原子,根据其原子序数获取对应特征向量
  • 堆叠所有原子的特征向量形成特征矩阵
  • 转换为PyTorch张量供深度学习模型使用

所以他的关键就是这一步

1
atom_fea = np.vstack([ari.get_atom_fea(crystal[i].specie.number)   for i in range(len(crystal))])

于是我就单独把它拎出来看一下

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
for i in range(len(crystal)):
print(crystal[i].specie.number)
'''
22
22
22
22
22
22
22
22
8
8
8
8
8
8
8
8
8
8
8
8
'''

正好是原子序数,这对我以后话题也很有启发

这样子画图就好很多了,而且能检索出元素了。

回归正题,在这里 crystal 变量是可以访问多个参数的,他是源自

1
from pymatgen.core.structure import Structure 

我们就来看这个 Structure 到底是什么

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
class Structure(IStructure, collections.abc.MutableSequence):
"""Mutable version of structure."""

__hash__ = None # type: ignore[assignment]

def __init__(
self,
lattice: ArrayLike | Lattice,
species: Sequence[CompositionLike],
coords: Sequence[ArrayLike] | ArrayLike,
charge: float | None = None,
validate_proximity: bool = False,
to_unit_cell: bool = False,
coords_are_cartesian: bool = False,
site_properties: dict | None = None,
labels: Sequence[str | None] | None = None,
properties: dict | None = None,
) -> None:
"""Create a periodic structure.

Args:
lattice: The lattice, either as a pymatgen.core.Lattice or
simply as any 2D array. Each row should correspond to a lattice
vector. e.g. [[10,0,0], [20,10,0], [0,0,30]] specifies a
lattice with lattice vectors [10,0,0], [20,10,0] and [0,0,30].
species: List of species on each site. Can take in flexible input,
including:

i. A sequence of element / species specified either as string
symbols, e.g. ["Li", "Fe2+", "P", ...] or atomic numbers,
e.g. (3, 56, ...) or actual Element or Species objects.

ii. List of dict of elements/species and occupancies, e.g.
[{"Fe" : 0.5, "Mn":0.5}, ...]. This allows the setup of
disordered structures.
coords (Nx3 array): list of fractional/cartesian coordinates of
each species.
charge (float): overall charge of the structure. Defaults to behavior
in SiteCollection where total charge is the sum of the oxidation
states.
validate_proximity (bool): Whether to check if there are sites
that are less than 0.01 Ang apart. Defaults to False.
to_unit_cell (bool): Whether to map all sites into the unit cell,
i.e., fractional coords between 0 and 1. Defaults to False.
coords_are_cartesian (bool): Set to True if you are providing
coordinates in Cartesian coordinates. Defaults to False.
site_properties (dict): Properties associated with the sites as a
dict of sequences, e.g. {"magmom":[5,5,5,5]}. The sequences
have to be the same length as the atomic species and
fractional_coords. Defaults to None for no properties.
labels (list[str]): Labels associated with the sites as a
list of strings, e.g. ['Li1', 'Li2']. Must have the same
length as the species and fractional coords. Defaults to
None for no labels.
properties (dict): Properties associated with the whole structure.
Will be serialized when writing the structure to JSON or YAML but is
lost when converting to other formats.
"""

参数:

  1. lattice: 晶格参数,可接受pymatgen.core.Lattice对象或任意二维数组。每行对应一个晶格向量,例如[[10,0,0], [20,10,0], [0,0,30]]表示晶格向量[10,0,0]、[20,10,0]和[0,0,30]。
  2. species: 各点位元素的种类列表,支持多种输入形式:
    1. 字符串符号(如[“Li”, “Fe2+”, “P”,…])、原子序数(如(3, 56,…))或Element/Species对象组成的序
    2. 由元素/物种及占位占比构成的字典列表(如[{“Fe”:0.5, “Mn”:0.5},…]),用于构建无序结构(disordered structures)coords (Nx3数组): 各元素的分数坐标/笛卡尔坐标列表
  3. charge (float): 结构整体电荷量。默认遵循SiteCollection行为,即总电荷为各点位氧化态之和
  4. validate_proximity (bool): 是否检查存在间距小于0.01埃的原子点位,默认为False
  5. to_unit_cell (bool): 是否将所有点位映射到单位晶胞内(即分数坐标强制转换到0-1区间),默认为False
  6. coords_are_cartesian (bool): True表示提供笛卡尔坐标,False表示分数坐标,默认为False
  7. site_properties (dict): 原子点位属性字典,格式为{“属性名”:[值1,值2,…]},序列长度需与species和coords一致,默认为None
  8. labels (list[str]): 原子点位标签列表(如[‘Li1’, ‘Li2’]),长度需与species和coords一致,默认为None
  9. properties (dict): 结构整体属性字典,序列化为JSON/YAML时保留,但转换为其他格式时会丢失

注:

  • 保留pymatgen.core.Lattice、Element、Species等专业对象名不翻译
  • 晶体学术语如”disordered structures”译为”无序结构”以匹配领域习惯
  • 坐标类型(笛卡尔/分数)的转换逻辑通过参数层次清晰表达
  • 默认值说明保持与API文档一致的数值精度(如0.01埃)
  • 数据结构格式(Nx3数组、字典、列表)保持技术文档通用表述
  • 参数间的数据一致性要求(如长度匹配)使用强调句式突显

Structure

以下是针对 pymatgen.core.structure.Structure 对象的详细访问指南,涵盖了其核心属性和方法:

1. 基本结构信息

1
2
3
4
5
6
7
8
9
10
11
# 晶胞参数
lattice = crystal.lattice
print(lattice.a, lattice.b, lattice.c) # 晶格常数 (Å)
print(lattice.alpha, lattice.beta, lattice.gamma) # 晶格角度 (°)
print(lattice.matrix) # 3x3 晶格向量矩阵

# 晶体基本信息
print(crystal.formula) # 化学式 (如 "Ti2 O3")
print(len(crystal)) # 原子总数
print(crystal.volume) # 晶胞体积 (ų)
print(crystal.density) # 密度 (g/cm³)

2. 原子位点访问

1
2
3
4
5
6
7
8
9
10
# 访问单个原子位点
site = crystal[0] # 第一个原子位点

# 原子位点属性
print(site.specie) # 元素对象 (如 Element("Ti"))
print(site.specie.number) # 原子序数 (22)
print(site.specie.symbol) # 元素符号 ("Ti")
print(site.coords) # 笛卡尔坐标 [x, y, z] (Å)
print(site.frac_coords) # 分数坐标 [a, b, c]
print(site.properties) # 附加属性字典

3. 化学组成分析

1
2
3
4
5
6
7
8
9
# 元素组成
print(crystal.composition) # Composition对象 (如 Composition("Ti2O3"))
print(crystal.composition.reduced_formula) # 最简化学式 ("Ti2O3")
print(crystal.composition.num_atoms) # 总原子数
print(crystal.composition["Ti"]) # Ti原子数量 (2)

# 元素分布
print(crystal.atomic_numbers) # 所有原子序数列表 [22, 22, 8, 8, 8]
print(crystal.species) # 所有元素对象列表

4. 对称性分析 (需要spglib)

1
2
3
4
5
6
7
8
# 空间群分析
print(crystal.get_space_group_info()) # (空间群符号, 国际表编号)

# 对称操作
print(crystal.get_symmetry_operations())

# 等效原子
print(crystal.find_equivalent_sites())

5. 结构操作

1
2
3
4
5
6
7
# 结构转换
crystal.to(fmt="poscar") # 转为VASP POSCAR格式字符串
crystal.to(filename="output.cif") # 保存为CIF文件

# 结构操作
crystal.make_supercell([2,2,2]) # 创建超胞
crystal.add_site_property("charge", [0.5]*len(crystal)) # 添加自定义属性

6. 物理性质计算

1
2
3
4
5
6
7
8
# 距离/角度计算
print(crystal.get_distance(0, 1)) # 原子0和1之间的距离
print(crystal.get_angle(0, 1, 2)) # 原子0-1-2之间的角度

# 邻居分析
neighbors = crystal.get_neighbors(site, r=3.0) # 3Å内的邻居
for neighbor, distance, index in neighbors:
print(f"原子{index}({neighbor.specie}) 距离: {distance:.3f}Å")

7. 高级特征

1
2
3
4
5
6
7
8
9
# 倒易晶格
print(crystal.lattice.reciprocal_lattice)

# 磁性分析 (若有自旋信息)
print(crystal.is_ordered) # 是否所有位点有序
print(crystal.spin_states) # 自旋状态

# 缺陷分析
print(crystal.frac_coords % 1) # 所有原子分数坐标

示例应用:获取所有原子坐标

1
2
3
4
5
6
7
8
9
10
11
# 获取所有原子的元素符号和笛卡尔坐标
atom_data = []
for i, site in enumerate(crystal):
atom_data.append({
"id": i,
"element": site.specie.symbol,
"atomic_number": site.specie.number,
"cartesian_coords": site.coords.tolist(),
"fractional_coords": site.frac_coords.tolist()
})
print(json.dumps(atom_data, indent=2))

熟悉了以上之后,便可以获取其原子磁矩了

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
for i in range(len(crystal)):
print(crystal[i].specie.symbol)
print(crystal.site_properties['magmom'][i])
'''
Ti
0.5139999999963609
Ti
-0.11199999999920497
Ti
0.9059999999935977
Ti
0.9059999999935977
Ti
0.5139999999963609
Ti
-0.11199999999920497
Ti
0.5139999999963609
Ti
0.5139999999963609
O
-0.006999999999953593
O
-0.006999999999953593
O
-0.006999999999953593
O
-0.006999999999953593
O
...
O
-0.006999999999953593
O
-0.006999999999953593
'''

现在就要把他加入进编码后的原子特征中,也就是 atom_fea。

首先把磁矩加入单个原子的特征中

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
atom_fea_list = []
for i in range(len(crystal)):
base_fea = ari.get_atom_fea(crystal[i].specie.number)
magmom = crystal.site_properties['magmom'][i] # 添加磁矩
magmom = float(magmom) # 确保磁矩是浮点数类型
# 将磁矩作为新特征追加
extended_fea = np.append(base_fea, magmom)

atom_fea_list.append(extended_fea)
print(extended_fea)
'''
[0. 0. 0. 0. 1. 0. 0. 0. 0. 0. 0. 0.
0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 1. 0.
0. 0. 0. 0. 1. 0. 0. 0. 0. 0. 0. 0.
0. 0. 0. 0. 1. 0. 0. 0. 0. 0. 0. 0.
0. 1. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0.
0. 1. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0.
1. 0. 0. 0. 0. 0. 0. 0. 1. 0. 0. 0.
0. 1. 0. 0. 0. 0. 0. 0. 0.514]
'''

再统一封装

1
2
3
4
5
6
7
8
9
10
11
12
13
14
atom_fea = np.vstack(atom_fea_list)
print('这是扩展后的原子特征')
atom_fea = torch.Tensor(atom_fea)
print(atom_fea)

'''
tensor([[ 0.0000, 0.0000, 0.0000, ..., 0.0000, 0.0000, 0.5140],
[ 0.0000, 0.0000, 0.0000, ..., 0.0000, 0.0000, -0.1120],
[ 0.0000, 0.0000, 0.0000, ..., 0.0000, 0.0000, 0.9060],
...,
[ 0.0000, 0.0000, 0.0000, ..., 0.0000, 0.0000, -0.0070],
[ 0.0000, 0.0000, 0.0000, ..., 0.0000, 0.0000, -0.0070],
[ 0.0000, 0.0000, 0.0000, ..., 0.0000, 0.0000, -0.0070]])
'''

关于磁矩特征的写入
http://baikelwang.github.io/2025/07/09/关于磁矩特征的写入/
作者
北海
发布于
2025年7月9日
许可协议