CN114818412A - 一种人脑受电磁辐射的高精度计算方法 - Google Patents

一种人脑受电磁辐射的高精度计算方法 Download PDF

Info

Publication number
CN114818412A
CN114818412A CN202210272635.1A CN202210272635A CN114818412A CN 114818412 A CN114818412 A CN 114818412A CN 202210272635 A CN202210272635 A CN 202210272635A CN 114818412 A CN114818412 A CN 114818412A
Authority
CN
China
Prior art keywords
absorption layer
calculation
magnetic field
electromagnetic
node
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Granted
Application number
CN202210272635.1A
Other languages
English (en)
Other versions
CN114818412B (zh
Inventor
陈广志
王新松
李尧尧
杜万里
苏东林
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Beihang University
Original Assignee
Beihang University
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Beihang University filed Critical Beihang University
Priority to CN202210272635.1A priority Critical patent/CN114818412B/zh
Publication of CN114818412A publication Critical patent/CN114818412A/zh
Application granted granted Critical
Publication of CN114818412B publication Critical patent/CN114818412B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F30/00—Computer-aided design [CAD]
    • G06F30/20—Design optimisation, verification or simulation
    • G06F30/23—Design optimisation, verification or simulation using finite element methods [FEM] or finite difference methods [FDM]
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F17/00—Digital computing or data processing equipment or methods, specially adapted for specific functions
    • G06F17/10—Complex mathematical operations
    • G06F17/11—Complex mathematical operations for solving equations, e.g. nonlinear equations, general mathematical optimization problems
    • G06F17/12—Simultaneous equations, e.g. systems of linear equations

Landscapes

  • Engineering & Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Mathematical Physics (AREA)
  • Theoretical Computer Science (AREA)
  • Mathematical Analysis (AREA)
  • Pure & Applied Mathematics (AREA)
  • Mathematical Optimization (AREA)
  • Computational Mathematics (AREA)
  • General Engineering & Computer Science (AREA)
  • Data Mining & Analysis (AREA)
  • Geometry (AREA)
  • Evolutionary Computation (AREA)
  • Computer Hardware Design (AREA)
  • Operations Research (AREA)
  • Algebra (AREA)
  • Databases & Information Systems (AREA)
  • Software Systems (AREA)
  • Measuring Magnetic Variables (AREA)
  • Management, Administration, Business Operations System, And Electronic Commerce (AREA)

Abstract

本发明公开了一种人脑受电磁辐射的高精度计算方法，包括以下步骤：S1.对于人脑模型进行定义:主要考虑人脑模型的尺寸，模拟脑模型内部的生物组织结构，并根据脑模型所处的电磁环境特点，赋予生物组织对应的电磁参数值；S2.人脑模型定义完成之后，使用FCC网格对于目标模型进行网格剖分：剖分过程中首先依据模型特点设置网格大小，然后确定每一个电磁场节点的位置；S3.使用FCC网格对于模型进行剖分后，对于场值迭代过程中参数进行初始化;S4.对吸收层内的参数初始化；S5.设置完吸收层的参数之后，设定时间步长
Figure DEST_PATH_IMAGE001
进行电磁场值的迭代计算。本发明对人脑暴露在电磁环境中可能存在的风险进行评估分析，并确保了人脑电磁辐射能量吸收的计算精度。

Description

一种人脑受电磁辐射的高精度计算方法
技术领域
本发明涉及电磁辐射，特别是涉及人脑受电磁辐射的高精度计算方法。
背景技术
针对日常生活中广泛存在的射频电磁环境，如使用手机上网、电话等情境将人体组织近距离长时间暴露在电磁辐射下，对人体健康造成了潜在的威胁。准确评估生物组织的电磁能量吸收对于衡量生物的电磁安全具有重要的意义。
时域有限差分(FDTD)方法利用差分方程处理麦克斯韦方程组中旋度方程，从时域角度出发解决电磁问题。因此FDTD方法在针对色散介质，如生物肌体组织、等离子体和雷达吸波材料时具有天然优势。传统FDTD方法计算精度主要受限于数值色散误差。基于面中心网格(FCC)的FCC-FDTD方法采用类似化学上SiC的晶胞结构，从基础的元胞角度出发，降低了传统FDTD方法的数值色散误差。然而为了使用FCC-FDTD方法评估脑模型受电磁辐射，需要开发一种高效率、高精度的吸收边界，从而达到在有限的计算内存中，使用FCC-FDTD方法分析电磁辐射对于人脑的潜在风险。
发明内容
本发明的目的在于克服现有技术的不足，提供一种人脑受电磁辐射的高精度计算方法，对人脑暴露在电磁环境中可能存在的风险进行评估分析，并确保了人脑电磁辐射能量吸收的计算精度。
本发明的目的是通过以下技术方案来实现的：一种人脑受电磁辐射的高精度计算方法，包括以下步骤：
S1.对于人脑模型进行定义:主要考虑人脑模型的尺寸，模拟脑模型内部的生物组织结构，并根据脑模型所处的电磁环境特点，赋予生物组织对应的电磁参数值；
进一步地，所述步骤S1中，依据电磁环境特点中的频率值，赋予生物组织对应的电磁参数值，所述电磁参数值包括介电常数值和磁导率值。
所述步骤S1中赋予生物组织对应的电磁参数值的过程包括：
S101.预先建立基于生物组织的频率与电磁参数对应表，表中包含多个频率，以及每一个频率对应的介电常数值和磁导率值；
S102.根据脑模型所处电磁环境的频率，从对应表中查找相应的介电常数值和磁导率值赋予生物组织。
S2.人脑模型定义完成之后，使用FCC网格对于目标模型进行网格剖分：剖分过程中首先依据模型特点设置网格大小，然后确定每一个电磁场节点的位置；
进一步地，所述步骤S2中，电磁场节点的位置包括电场节点坐标和磁场节点坐标；
设网格大小为Δx、Δy、Δz，FCC网格坐标为(isΔx，jsΔy，ksΔz)，其中Δi(i＝x，y，z)为FCC网格在i方向上的大小；is，js，ks是当前网格在x、y、z三个方向上的编号；
其中分布的四种电场节点坐标分别为：
E1＝(isΔx，jsΔy，ksΔz)
E2＝((is+0.5)Δx，(js+0.5)Δy，ksΔz)
E3＝((is+0.5)Δx，jsΔy，(ks+0.5)Δz)
E4＝(isΔx，(js+0.5)Δy，(ks+0.5)Δz)
该网格内分布的四种磁场节点坐标分别为：
Figure BDA0003554316030000021
S3.使用FCC网格对于模型进行剖分后，对于场值迭代过程中参数进行初始化；
进一步地，所述步骤S3中的初始化内容包括：将整个计算区域的四种电场值和四种磁场值初始化为零；
并依据人脑模型所处环境的频率值，确定辐射源类型：
预先建立频率值与辐射源类型的对应表，表中包含多个频率，以及每一个频率对应的辐射源类型；每一类辐射源均具有已知的时域波形和激励时长；
根据人脑模型所处环境的频率值，从频率值与辐射源参数的对应表中查找相应的辐射源类型，并确定该类辐射源的时域波形和激励时长。
S4.在计算区域初始化之后，对吸收层内的参数初始化。
进一步地，所述步骤S4的关键在于确定吸收层内每个网格中每个电磁场节点到吸收层和计算区域分界面的距离；
所述步骤S4包括以下子步骤：
S401.初始化吸收层内磁场节点到吸收层和分界面的距离；四类磁场节点到吸收层与计算区域的分界面有两类距离，分别为ρm，1(k)＝(k-0.25)Δi，ρm，2(k)＝(k-0.75)Δi，其中ρm，1(k)、ρm，2(k)表示从磁场节点到吸收层和计算区域分界面的距离，k表示磁场节点到吸收层和计算区域分界面的网格数；
S402.初始化吸收层内电场节点到吸收层和分界面的距离；四类电场节点到吸收层与计算区域的分界面有两种距离，分别为ρe，1(k)＝kΔi，ρe，2(k)＝(k-0.5)Δi，ρe，1(k)、ρe，2(k)表示电场节点到吸收层和计算区域分界面的距离，其中k表示电场节点到吸收层和计算区域分界面的网格数。
S5.设置完吸收层的参数之后，设定时间步长Δt，进行电磁场值的迭代计算。
进一步地，所述步骤S5包括：
S501.设定时间步长Δt，n为当前时刻迭代步数，在(n+0.5)Δt时刻，更新计算区域和吸收层内的磁场值，对于第一类磁场节点H1，其迭代计算公式计算如下：
Figure BDA0003554316030000031
Figure BDA0003554316030000032
Figure BDA0003554316030000033
其中μ为磁导率，
Figure BDA0003554316030000034
和
Figure BDA0003554316030000035
和的计算方式取决于FCC网格中磁场节点H1周围环绕的电场节点分布方式；
四类磁场节点H1～H4的计算迭代公式仅存在角标的区别，对于H2～H4进行迭代时，按照H1的计算方式进行，并将下角标进行替换；
S502.在更新完计算区域和吸收层内的磁场值之后，对吸收层的磁场进行再次更新；当x方向设有吸收层时，需要对x方向吸收层内Hy，Hz进行更新计算；当y方向设有吸收层时，需要对y方向吸收层内Hx，Hz进行更新计算；当z方向设有吸收层时，需要对z方向吸收层内Hx，Hy进行更新计算，依据给定的磁场节点；
吸收层磁场节点H1～H4的计算迭代公式仅存在角标的区别，对于第一类磁场节点H1，其迭代计算公式计算如下：
Figure BDA0003554316030000036
Figure BDA0003554316030000037
Figure BDA0003554316030000038
κ为吸收层参数；对于吸收层内H2～H4的计算方式按照H1的计算方式进行，并在计算时将下角标进行替换，
Figure BDA0003554316030000039
和
Figure BDA00035543160300000310
是吸收层内针对第一类磁场节点H1迭代更新的过程量；吸收层内H2～H4迭代更新的过程量为
Figure BDA00035543160300000311
依据标注的磁场节点坐标，及磁场节点到吸收层和计算区域分别面，将
Figure BDA00035543160300000312
划分为两类，并由ρm，1和ρm，2决定；在x方向上的吸收层内
Figure BDA00035543160300000313
和
Figure BDA00035543160300000314
由ρm，1和磁场节点类型决定，
Figure BDA00035543160300000315
和
Figure BDA00035543160300000316
由ρm，2和磁场节点类型决定；y方向上的吸收层内
Figure BDA0003554316030000041
和
Figure BDA0003554316030000042
由ρm，1和磁场节点类型决定，
Figure BDA0003554316030000043
和
Figure BDA0003554316030000044
由ρm，2和磁场节点类型决定；z方向上的吸收层内
Figure BDA0003554316030000045
和
Figure BDA0003554316030000046
由ρm，1和磁场节点类型决定，
Figure BDA0003554316030000047
和
Figure BDA0003554316030000048
由ρm，2和磁场节点类型决定；
以ρm，1和第一类磁场节点决定的
Figure BDA0003554316030000049
的计算公式进行说明：
Figure BDA00035543160300000410
其中αm，1，σm，1和κm，1是吸收层参数，并由ρm，1决定，计算方法如下
Figure BDA00035543160300000411
Figure BDA00035543160300000412
Figure BDA00035543160300000413
其中σmax，αmax，κmax和ncpml是吸收层固定常数，δ是吸收层网格数乘此网格步长；
由于计算方式相似，
Figure BDA00035543160300000414
和
Figure BDA00035543160300000415
相对于
Figure BDA00035543160300000416
仅需要修改下角标的方向，
Figure BDA00035543160300000417
仅需要对于下角标编号更换，就能够完成计算；
S503.处理完磁场之后，在(n+1)Δt时刻更新计算区域和吸收层内的电场值；
Figure BDA00035543160300000418
Figure BDA00035543160300000419
Figure BDA00035543160300000420
ε为介电常数，
Figure BDA00035543160300000421
和
Figure BDA00035543160300000422
的计算方式取决于FCC网格中，当前迭代的电场节点周围环绕的磁场节点分布方式；
S504.在更新完计算区域和吸收层内的电场值之后，对吸收层的电场进行再次更新；针对两种距离的磁场节点类型；当x方向设有吸收层时，需要对吸收层内Ey，Ez进行更新计算；当y方向设有吸收层时，需要对吸收层内Ex，Ez进行更新计算；当z方向设有吸收层时，需要对吸收层内Ex，Ey进行更新计算；对于第一类电场节点E1，其迭代计算方式如下：
Figure BDA00035543160300000423
Figure BDA00035543160300000424
Figure BDA00035543160300000425
对于吸收层内E2～E4的计算方式按照E1的计算方式，并在下角标进行替换，
Figure BDA0003554316030000051
和
Figure BDA0003554316030000052
是吸收层内针对第一类磁场节点E1迭代更新的过程量；吸收层内E2～E4迭代更新的过程量为
Figure BDA0003554316030000053
依据标注的电场节点坐标，及电场节点到吸收层和计算区域分界面，将
Figure BDA0003554316030000054
分为两类，并由ρe，1和ρe，2决定；在x方向上的吸收层内
Figure BDA0003554316030000055
和
Figure BDA0003554316030000056
由ρe，1和电场节点类型决定，
Figure BDA0003554316030000057
和
Figure BDA0003554316030000058
由ρe，2和电场节点类型决定；y方向上的吸收层内
Figure BDA0003554316030000059
和
Figure BDA00035543160300000510
由ρe，1和电场节点类型决定，
Figure BDA00035543160300000511
和
Figure BDA00035543160300000512
由ρe，2和电场节点类型决定；z方向上的吸收层内
Figure BDA00035543160300000513
和
Figure BDA00035543160300000514
由ρe，1和电场节点类型决定，
Figure BDA00035543160300000515
和
Figure BDA00035543160300000516
由ρe，2和电场节点类型决定；
以ρe，1和第一类电场节点决定的
Figure BDA00035543160300000517
的计算公式进行说明：
Figure BDA00035543160300000518
其中αe，1，σe，1和κe，1是吸收层参数，并由ρe，1决定，计算方法如下
Figure BDA00035543160300000519
Figure BDA00035543160300000520
Figure BDA00035543160300000521
由于计算方式相似，
Figure BDA00035543160300000522
和
Figure BDA00035543160300000523
相对于
Figure BDA00035543160300000524
仅需要修改下角标的方向，
Figure BDA00035543160300000525
仅需要对于下角标编号更换，就能够完成计算；
S505.将时间步由n时刻变成n+1时刻，并进行当前时间步迭代结果输出，输出当前时刻记录点处的电磁场值，即人脑受到的电磁辐射；
判断是否达到设置的迭代步数最大值，若否，返回步骤S501，进行新一轮的电磁场值计算，若是，计算停止。
本发明的有益效果是：本发明使用基于FCC-CPML吸收边界的FCC-FDTD方法分析人脑模型受电磁辐射，从FCC-FDTD方法出发，推导了适用于FCC-FDTD方法的FCC-CPML边界。FCC-CPML边界可以确保FCC-FDTD方法使用有限的计算空间，在长时间的人脑模型受电磁辐射的计算分析中，电磁波无反射透过边界，确保了人脑电磁辐射能量吸收的计算精度，也就是说，本发明能够精准计算人脑模型受电磁辐射时，人脑模型内部生物组织吸收的电磁能量。
附图说明
图1为本发明的方法流程图；
图2实施例中人脑模型及网格剖分示意图；
图3为实施例中FDTD方法计算人脑模型吸收的电磁能量示意图；
图4为实施例中基于FCC-CPML吸收边界的FCC-FDTD方法计算人脑模型吸收的电磁能量。
具体实施方式
下面结合附图进一步详细描述本发明的技术方案，但本发明的保护范围不局限于以下所述。
如图1所示，一种人脑受电磁辐射的高精度计算方法，包括以下步骤：
S1.对于人脑模型进行定义：主要考虑人脑模型的尺寸，模拟脑模型内部的生物组织结构，并根据脑模型所处的电磁环境特点，赋予生物组织对应的电磁参数值；
进一步地，所述步骤S1中，依据电磁环境特点中的频率值，赋予生物组织对应的电磁参数值，所述电磁参数值包括介电常数值和磁导率值。
所述步骤S1中赋予生物组织对应的电磁参数值的过程包括：
S101.预先建立基于生物组织的频率与电磁参数对应表，表中包含多个频率，以及每一个频率对应的介电常数值和磁导率值；
S102.根据脑模型所处电磁环境的频率，从对应表中查找相应的介电常数值和磁导率值赋予生物组织。
S2.人脑模型定义完成之后，使用FCC网格对于目标模型进行网格剖分：剖分过程中首先依据模型特点设置网格大小，然后确定每一个电磁场节点的位置；
进一步地，所述步骤S2中，电磁场节点的位置包括电场节点坐标和磁场节点坐标；
设网格大小为Δx、Δy、Δz，FCC网格坐标为(isΔx，jsΔy，ksΔz)，其中Δi(i＝x，y，z)为FCC网格在i方向上的大小；is，js，ks是当前网格在x、y、z三个方向上的编号；
其中分布的四种电场节点坐标分别为：
E1＝(isΔx，jsΔy，ksΔz)
E2＝((is+0.5)Δx，(js+0.5)Δy，ksΔz)
E3＝((iS+0.5)Δx，jsΔy，(ks+0.5)Δz)
E4＝(isΔx，(js+0.5)Δy，(ks+0.5)Δz)
该网格内分布的四种磁场节点坐标分别为：
Figure BDA0003554316030000071
S3.使用FCC网格对于模型进行剖分后，对于场值迭代过程中参数进行初始化；
进一步地，所述步骤S3中的初始化内容包括：将整个计算区域的四种电场值和四种磁场值初始化为零；
并依据人脑模型所处环境的频率值，确定辐射源类型：
预先建立频率值与辐射源类型的对应表，表中包含多个频率，以及每一个频率对应的辐射源类型；每一类辐射源均具有已知的时域波形和激励时长；
根据人脑模型所处环境的频率值，从频率值与辐射源参数的对应表中查找相应的辐射源类型，并确定该类辐射源的时域波形和激励时长。
S4.在计算区域初始化之后，对吸收层内的参数初始化。
进一步地，所述步骤S4的关键在于确定吸收层内每个网格中每个电磁场节点到吸收层和计算区域分界面的距离；
所述步骤S4包括以下子步骤：
S401.初始化吸收层内磁场节点到吸收层和分界面的距离；四类磁场节点到吸收层与计算区域的分界面有两类距离，分别为ρm，1(k)＝(k-0.25)Δi，ρm，2(k)＝(k-0.75)Δi，其中ρm，1(k)、ρm，2(k)表示从磁场节点到吸收层和计算区域分界面的距离，k表示磁场节点到吸收层和计算区域分界面的网格数；
S402.初始化吸收层内电场节点到吸收层和分界面的距离；四类电场节点到吸收层与计算区域的分界面有两种距离，分别为ρe，1(k)＝kΔi，ρe，2(k)＝(k-0.5)Δi，ρe，1(k)、ρe，2(k)表示电场节点到吸收层和计算区域分界面的距离，其中k表示电场节点到吸收层和计算区域分界面的网格数。
S5.设置完吸收层的参数之后，设定时间步长Δt，进行电磁场值的迭代计算。
进一步地，所述步骤S5包括：
S501.设定时间步长Δt，n为当前时刻迭代步数，在(n+0.5)Δt时刻，更新计算区域和吸收层内的磁场值，对于第一类磁场节点H1，其迭代计算公式计算如下：
Figure BDA0003554316030000072
Figure BDA0003554316030000073
Figure BDA0003554316030000081
其中μ为磁导率，
Figure BDA0003554316030000082
和
Figure BDA0003554316030000083
和的计算方式取决于FCC网格中磁场节点H1周围环绕的电场节点分布方式；
四类磁场节点H1～H4的计算迭代公式仅存在角标的区别，对于H2～H4进行迭代时，按照H1的计算方式进行，并将下角标进行替换；
S502.在更新完计算区域和吸收层内的磁场值之后，对吸收层的磁场进行再次更新；当x方向设有吸收层时，需要对x方向吸收层内Hy，Hz进行更新计算；当y方向设有吸收层时，需要对y方向吸收层内Hx，Hz进行更新计算；当z方向设有吸收层时，需要对z方向吸收层内Hx，Hy进行更新计算，依据给定的磁场节点；
吸收层磁场节点H1～H4的计算迭代公式仅存在角标的区别，对于第一类磁场节点H1，其迭代计算公式计算如下：
Figure BDA0003554316030000084
Figure BDA0003554316030000085
Figure BDA0003554316030000086
κ为吸收层参数；对于吸收层内H2～H4的计算方式按照H1的计算方式进行，并在计算时将下角标进行替换，
Figure BDA0003554316030000087
和
Figure BDA0003554316030000088
是吸收层内针对第一类磁场节点H1迭代更新的过程量；吸收层内H2～H4迭代更新的过程量为
Figure BDA0003554316030000089
依据标注的磁场节点坐标，及磁场节点到吸收层和计算区域分别面，将
Figure BDA00035543160300000810
划分为两类，并由ρm，1和ρm，2决定；在x方向上的吸收层内
Figure BDA00035543160300000811
和
Figure BDA00035543160300000812
由ρm，1和磁场节点类型决定，
Figure BDA00035543160300000813
和
Figure BDA00035543160300000814
由ρm，2和磁场节点类型决定；y方向上的吸收层内
Figure BDA00035543160300000815
和
Figure BDA00035543160300000816
由ρm，1和磁场节点类型决定，
Figure BDA00035543160300000817
和
Figure BDA00035543160300000818
由ρm，2和磁场节点类型决定；z方向上的吸收层内
Figure BDA00035543160300000819
和
Figure BDA00035543160300000820
由ρm，1和磁场节点类型决定，
Figure BDA00035543160300000821
和
Figure BDA00035543160300000822
由ρm，2和磁场节点类型决定；
以ρm，1和第一类磁场节点决定的
Figure BDA00035543160300000823
的计算公式进行说明：
Figure BDA00035543160300000824
其中αm，1，σm，1和κm，1是吸收层参数，并由ρm，1决定，计算方法如下
Figure BDA0003554316030000091
Figure BDA0003554316030000092
Figure BDA0003554316030000093
其中σmax，αmax，κmax和ncpml是吸收层固定常数，δ是吸收层网格数乘此网格步长；
由于计算方式相似，
Figure BDA0003554316030000094
和
Figure BDA0003554316030000095
相对于
Figure BDA0003554316030000096
仅需要修改下角标的方向，
Figure BDA0003554316030000097
仅需要对于下角标编号更换，就能够完成计算；
S503.处理完磁场之后，在(n+1)Δt时刻更新计算区域和吸收层内的电场值；
Figure BDA0003554316030000098
Figure BDA0003554316030000099
Figure BDA00035543160300000910
ε为介电常数，
Figure BDA00035543160300000911
和
Figure BDA00035543160300000912
的计算方式取决于FCC网格中，当前迭代的电场节点周围环绕的磁场节点分布方式；
S504.在更新完计算区域和吸收层内的电场值之后，对吸收层的电场进行再次更新；针对两种距离的磁场节点类型；当x方向设有吸收层时，需要对吸收层内Ey，Ez进行更新计算；当y方向设有吸收层时，需要对吸收层内Ex，Ez进行更新计算；当z方向设有吸收层时，需要对吸收层内Ex，Ey进行更新计算；对于第一类电场节点E1，其迭代计算方式如下：
Figure BDA00035543160300000913
Figure BDA00035543160300000914
Figure BDA00035543160300000915
对于吸收层内E2～E4的计算方式按照E1的计算方式，并在下角标进行替换，
Figure BDA00035543160300000916
和
Figure BDA00035543160300000917
是吸收层内针对第一类磁场节点E1迭代更新的过程量；吸收层内E2～E4迭代更新的过程量为
Figure BDA00035543160300000918
依据标注的电场节点坐标，及电场节点到吸收层和计算区域分界面，将
Figure BDA00035543160300000919
分为两类，并由ρe，1和ρe，2决定；在x方向上的吸收层内
Figure BDA00035543160300000920
和
Figure BDA00035543160300000921
由ρe，1和电场节点类型决定，
Figure BDA00035543160300000922
和
Figure BDA00035543160300000923
由ρe，2和电场节点类型决定；y方向上的吸收层内
Figure BDA00035543160300000924
和
Figure BDA00035543160300000925
由ρe，1和电场节点类型决定，
Figure BDA00035543160300000926
和
Figure BDA0003554316030000101
由ρe，2和电场节点类型决定；z方向上的吸收层内
Figure BDA0003554316030000102
和
Figure BDA0003554316030000103
由ρe，1和电场节点类型决定，
Figure BDA0003554316030000104
和
Figure BDA0003554316030000105
由ρe，2和电场节点类型决定；
以ρe，1和第一类电场节点决定的
Figure BDA0003554316030000106
的计算公式进行说明：
Figure BDA0003554316030000107
其中αe，1，σe，1和κe，1是吸收层参数，并由ρe，1决定，计算方法如下
Figure BDA0003554316030000108
Figure BDA0003554316030000109
Figure BDA00035543160300001010
由于计算方式相似，
Figure BDA00035543160300001011
和
Figure BDA00035543160300001012
相对于
Figure BDA00035543160300001013
仅需要修改下角标的方向，
Figure BDA00035543160300001014
仅需要对于下角标编号更换，就能够完成计算；
S505.将时间步由n时刻变成n+1时刻，并进行当前时间步迭代结果输出，输出当前时刻记录点处的电磁场值，即人脑受到的电磁辐射；
判断是否达到设置的迭代步数最大值，若否，返回步骤S501，进行新一轮的电磁场值计算，若是，计算停止。
下面结合具体实施例对本申请作进一步说明：
图2表示的是一个考虑了四层介质的人脑模型，分析当手机使用900MHz进行通信时，人脑模型在一个时域周期内吸收的电磁能量。该人脑模型由四种生物组织组成，分别是皮肤(厚度：4cm)，脂肪(厚度：4cm)，骨头(厚度：8cm)和一种人脑内部生物组织进行平均后的生物组织(厚度：92cm)。考虑使用2cm×2cm×2cm的FCC网格大小对人脑模型进行剖分。利用本申请提供的一种人脑受电磁辐射高精度计算方法的具体步骤如下：
步骤一：对于目标进行定义、剖分和FCC网格的初始化。
步骤二：确定手机使用900MHz通信时的辐射源形式。为了更好的模拟手机使用900MHz进行通信时，本实施例使用如下辐射源代替手机
Figure BDA00035543160300001015
式中τ＝t/(4π/ω)。该辐射源先经历两个渐变周期，到达稳态，并在稳态运行四个周期后计算最后一个周期内脑模型吸收的电磁能量。
步骤三：基于FCC网格，对麦克斯韦方程组进行时域有限差分变换处理，得到空间中电磁场的分布。
步骤四：对于获得的计算周期内的电磁场值进行后处理。在获得了一个辐射源辐射周期内的电磁场值之后，我们利用如下公式：
Figure BDA0003554316030000111
获得人脑模型在一个辐射周期内吸收的电磁能量。图3是使用传统FDTD方法获得的在一个辐射周期内的人脑模型吸收的电磁能量，脑模型吸收的电磁能量最大值为Pmax＝0.2453W/m3。图4是使用基于FCC-CPML吸收边界的FCC-FDTD方法获得的在一个辐射周期内人脑模型吸收生物电磁能量，脑模型吸收的电磁能量最大值为Pmax＝0.2024W/m3。使用幂级数计算出脑模型吸收的电磁能量最大值的精确结果为Pmax＝0.1527W/m3。使用本专利提出的方法相对于传统方法提升的计算精度为46.33％。
上述说明示出并描述了本发明的一个优选实施例，但如前所述，应当理解本发明并非局限于本文所披露的形式，不应看作是对其他实施例的排除，而可用于各种其他组合、修改和环境，并能够在本文所述发明构想范围内，通过上述教导或相关领域的技术或知识进行改动。而本领域人员所进行的改动和变化不脱离本发明的精神和范围，则都应在本发明所附权利要求的保护范围内。

Claims (7)

1.一种人脑受电磁辐射的高精度计算方法，其特征在于：包括以下步骤：
S1.对于人脑模型进行定义：主要考虑人脑模型的尺寸，模拟脑模型内部的生物组织结构，并根据脑模型所处的电磁环境特点，赋予生物组织对应的电磁参数值；
S2.人脑模型定义完成之后，使用FCC网格对于目标模型进行网格剖分：剖分过程中首先依据模型特点设置网格大小，然后确定每一个电磁场节点的位置；
S3.使用FCC网格对于模型进行剖分后，对于场值迭代过程中参数进行初始化；
S4.在计算区域初始化之后，对吸收层内的参数初始化；
S5.设置完吸收层的参数之后，设定时间步长Δt，进行电磁场值的迭代计算。
2.根据权利要求1所述的一种人脑受电磁辐射的高精度计算方法，其特征在于：所述步骤S1中，依据电磁环境特点中的频率值，赋予生物组织对应的电磁参数值，所述电磁参数值包括介电常数值和磁导率值。
3.根据权利要求2所述的一种人脑受电磁辐射的高精度计算方法，其特征在于：所述步骤S1中赋予生物组织对应的电磁参数值的过程包括：
S101.预先建立基于生物组织的频率与电磁参数对应表，表中包含多个频率，以及每一个频率对应的介电常数值和磁导率值；
S102.根据脑模型所处电磁环境的频率，从对应表中查找相应的介电常数值和磁导率值赋予生物组织。
4.根据权利要求1所述的一种人脑受电磁辐射的高精度计算方法，其特征在于：所述步骤S2中，电磁场节点的位置包括电场节点坐标和磁场节点坐标；
设网格大小为Δx、Δy、Δz，FCC网格坐标为(isΔx，jsΔy，ksΔz)，其中Δi(i＝x，y，z)为FCC网格在i方向上的大小；is，js，ks是当前网格在x、y、z三个方向上的编号；
其中分布的四种电场节点坐标分别为：
E1＝(isΔx，jsΔy，ksΔz)
E2＝((is+0.5)Δx，(js+0.5)Δy，ksΔz)
E3＝((is+0.5)Δx，jsΔy，(ks+0.5)Δz)
E4＝(isΔx，(js+0.5)Δy，(ks+0.5)Δz)
该网格内分布的四种磁场节点坐标分别为：
Figure FDA0003554316020000011
5.根据权利要求1所述的一种人脑受电磁辐射的高精度计算方法，其特征在于：所述步骤S3中的初始化内容包括：将整个计算区域的四种电场值和四种磁场值初始化为零；
并依据人脑模型所处环境的频率值，确定辐射源类型：
预先建立频率值与辐射源类型的对应表，表中包含多个频率，以及每一个频率对应的辐射源类型；每一类辐射源均具有已知的时域波形和激励时长；
根据人脑模型所处环境的频率值，从频率值与辐射源参数的对应表中查找相应的辐射源类型，并确定该类辐射源的时域波形和激励时长。
6.根据权利要求1所述的一种人脑受电磁辐射的高精度计算方法，其特征在于：所述步骤S4的关键在于确定吸收层内每个网格中每个电磁场节点到吸收层和计算区域分界面的距离；
所述步骤S4包括以下子步骤：
S401.初始化吸收层内磁场节点到吸收层和分界面的距离；四类磁场节点到吸收层与计算区域的分界面有两类距离，分别为ρm，1(k)＝(k-0.25)Δi，ρm，2(k)＝(k-0.75)Δi，其中ρm，1(k)、ρm，2(k)表示从磁场节点到吸收层和计算区域分界面的距离，k表示磁场节点到吸收层和计算区域分界面的网格数；
S402.初始化吸收层内电场节点到吸收层和分界面的距离；四类电场节点到吸收层与计算区域的分界面有两种距离，分别为ρe，1(k)＝kΔi，ρe，2(k)＝(k-0.5)Δi，ρe，1(k)、ρe，2(k)表示电场节点到吸收层和计算区域分界面的距离，其中k表示电场节点到吸收层和计算区域分界面的网格数。
7.根据权利要求1所述的一种人脑受电磁辐射的高精度计算方法，其特征在于：所述步骤S5包括：
S501.设定时间步长Δt，n为当前时刻迭代步数，在(n+0.5)Δt时刻，更新计算区域和吸收层内的磁场值，对于第一类磁场节点H1，其迭代计算公式计算如下：
Figure FDA0003554316020000021
Figure FDA0003554316020000022
Figure FDA0003554316020000023
其中μ为磁导率，
Figure FDA0003554316020000024
和
Figure FDA0003554316020000025
和的计算方式取决于FCC网格中磁场节点H1周围环绕的电场节点分布方式；
四类磁场节点H1～H4的计算迭代公式仅存在角标的区别，对于H2～H4进行迭代时，按照H1的计算方式进行，并将下角标进行替换；
S502.在更新完计算区域和吸收层内的磁场值之后，对吸收层的磁场进行再次更新；当x方向设有吸收层时，需要对x方向吸收层内Hy，Hz进行更新计算；当y方向设有吸收层时，需要对y方向吸收层内Hx，Hz进行更新计算；当z方向设有吸收层时，需要对z方向吸收层内Hx，Hy进行更新计算，依据给定的磁场节点；
吸收层磁场节点H1～H4的计算迭代公式仅存在角标的区别，对于第一类磁场节点H1，其迭代计算公式计算如下：
Figure FDA0003554316020000031
Figure FDA0003554316020000032
Figure FDA0003554316020000033
κ为吸收层参数；对于吸收层内H2～H4的计算方式按照H1的计算方式进行，并在计算时将下角标进行替换，
Figure FDA0003554316020000034
和
Figure FDA0003554316020000035
是吸收层内针对第一类磁场节点H1迭代更新的过程量；吸收层内H2～H4迭代更新的过程量为
Figure FDA0003554316020000036
依据标注的磁场节点坐标，及磁场节点到吸收层和计算区域分别面，将
Figure FDA0003554316020000037
划分为两类，并由ρm，1和ρm，2决定；在x方向上的吸收层内
Figure FDA0003554316020000038
和
Figure FDA0003554316020000039
由ρm，1和磁场节点类型决定，
Figure FDA00035543160200000310
和
Figure FDA00035543160200000311
由ρm，2和磁场节点类型决定；y方向上的吸收层内
Figure FDA00035543160200000312
和
Figure FDA00035543160200000313
由ρm，1和磁场节点类型决定，
Figure FDA00035543160200000314
和
Figure FDA00035543160200000315
由ρm，2和磁场节点类型决定；z方向上的吸收层内
Figure FDA00035543160200000316
和
Figure FDA00035543160200000317
由ρm，1和磁场节点类型决定，
Figure FDA00035543160200000318
和
Figure FDA00035543160200000319
由ρm，2和磁场节点类型决定；
以ρm，1和第一类磁场节点决定的
Figure FDA00035543160200000320
的计算公式进行说明：
Figure FDA00035543160200000321
其中αm，1，σm，1和κm，1是吸收层参数，并由ρm，1决定，计算方法如下
Figure FDA00035543160200000322
Figure FDA00035543160200000323
Figure FDA00035543160200000324
其中σmax，αmax，κmax和ncpml是吸收层固定常数，δ是吸收层网格数乘此网格步长；
由于计算方式相似，
Figure FDA00035543160200000325
和
Figure FDA00035543160200000326
相对于
Figure FDA00035543160200000327
仅需要修改下角标的方向，
Figure FDA00035543160200000328
仅需要对于下角标编号更换，就能够完成计算；
S503.处理完磁场之后，在(n+1)Δt时刻更新计算区域和吸收层内的电场值；
Figure FDA0003554316020000041
Figure FDA0003554316020000042
Figure FDA0003554316020000043
ε为介电常数，
Figure FDA0003554316020000044
和
Figure FDA0003554316020000045
的计算方式取决于FCC网格中，当前迭代的电场节点周围环绕的磁场节点分布方式；
S504.在更新完计算区域和吸收层内的电场值之后，对吸收层的电场进行再次更新；针对两种距离的磁场节点类型；当x方向设有吸收层时，需要对吸收层内Ey，Ez进行更新计算；当y方向设有吸收层时，需要对吸收层内Ex，Ez进行更新计算；当z方向设有吸收层时，需要对吸收层内Ex，Ey进行更新计算；对于第一类电场节点E1，其迭代计算方式如下：
Figure FDA0003554316020000046
Figure FDA0003554316020000047
Figure FDA0003554316020000048
对于吸收层内E2～E4的计算方式按照E1的计算方式，并在下角标进行替换，
Figure FDA0003554316020000049
和
Figure FDA00035543160200000410
是吸收层内针对第一类磁场节点E1迭代更新的过程量；吸收层内E2～E4迭代更新的过程量为
Figure FDA00035543160200000411
依据标注的电场节点坐标，及电场节点到吸收层和计算区域分界面，将
Figure FDA00035543160200000412
分为两类，并由ρe，1和ρe，2决定；在x方向上的吸收层内
Figure FDA00035543160200000413
和
Figure FDA00035543160200000414
由ρe，1和电场节点类型决定，
Figure FDA00035543160200000415
和
Figure FDA00035543160200000416
由ρe，2和电场节点类型决定；y方向上的吸收层内
Figure FDA00035543160200000417
和
Figure FDA00035543160200000418
由ρe，1和电场节点类型决定，
Figure FDA00035543160200000419
和
Figure FDA00035543160200000420
由ρe，2和电场节点类型决定；z方向上的吸收层内
Figure FDA00035543160200000421
和
Figure FDA00035543160200000422
由ρe，1和电场节点类型决定，
Figure FDA00035543160200000423
和
Figure FDA00035543160200000424
由ρe，2和电场节点类型决定；
以ρe，1和第一类电场节点决定的
Figure FDA00035543160200000425
的计算公式进行说明：
Figure FDA00035543160200000426
其中αe，1，σe，1和κe，1是吸收层参数，并由ρe，1决定，计算方法如下
Figure FDA0003554316020000051
Figure FDA0003554316020000052
Figure FDA0003554316020000053
由于计算方式相似，
Figure FDA0003554316020000054
和
Figure FDA0003554316020000055
相对于
Figure FDA0003554316020000056
仅需要修改下角标的方向，
Figure FDA0003554316020000057
仅需要对于下角标编号更换，就能够完成计算；
S505.将时间步由n时刻变成n+1时刻，并进行当前时间步迭代结果输出，输出当前时刻记录点处的电磁场值，即人脑受到的电磁辐射；
判断是否达到设置的迭代步数最大值，若否，返回步骤S501，进行新一轮的电磁场值计算，若是，计算停止。
CN202210272635.1A 2022-03-18 2022-03-18 一种人脑受电磁辐射的高精度计算方法 Active CN114818412B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202210272635.1A CN114818412B (zh) 2022-03-18 2022-03-18 一种人脑受电磁辐射的高精度计算方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202210272635.1A CN114818412B (zh) 2022-03-18 2022-03-18 一种人脑受电磁辐射的高精度计算方法

Publications (2)

Publication Number Publication Date
CN114818412A true CN114818412A (zh) 2022-07-29
CN114818412B CN114818412B (zh) 2022-09-27

Family

ID=82530811

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202210272635.1A Active CN114818412B (zh) 2022-03-18 2022-03-18 一种人脑受电磁辐射的高精度计算方法

Country Status (1)

Country Link
CN (1) CN114818412B (zh)

Cited By (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN117236110A (zh) * 2023-04-13 2023-12-15 安徽大学 基于电热多物理场随机fdtd的电磁辐射强度评估方法

Citations (8)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN1473268A (zh) * 2001-08-08 2004-02-04 株式会社Ntt都科摩 吸收功率测量装置
CN102013895A (zh) * 2010-12-13 2011-04-13 惠州市硕贝德通讯科技有限公司 一种解决小型终端手机天线的电磁兼容性方法
CN102437410A (zh) * 2011-09-19 2012-05-02 惠州硕贝德无线科技股份有限公司 提升手机天线与助听器电磁兼容性的方法
CN102510296A (zh) * 2011-11-09 2012-06-20 中兴通讯股份有限公司 移动终端及降低移动终端辐射的方法
CN108027412A (zh) * 2015-07-07 2018-05-11 Q生物公司 场不变定量磁共振特征标志
CN108136195A (zh) * 2015-10-26 2018-06-08 奥海视网膜科技有限公司 使用脉冲的能量源进行生物组织热疗的方法
CN110998971A (zh) * 2017-05-23 2020-04-10 撒加德有限公司 辐射屏蔽罩
WO2022031362A1 (en) * 2020-08-04 2022-02-10 Ojai Retinal Technology, Llc System and method for preventing or treating alzheimer's and other neurodegenerative diseases

Patent Citations (8)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN1473268A (zh) * 2001-08-08 2004-02-04 株式会社Ntt都科摩 吸收功率测量装置
CN102013895A (zh) * 2010-12-13 2011-04-13 惠州市硕贝德通讯科技有限公司 一种解决小型终端手机天线的电磁兼容性方法
CN102437410A (zh) * 2011-09-19 2012-05-02 惠州硕贝德无线科技股份有限公司 提升手机天线与助听器电磁兼容性的方法
CN102510296A (zh) * 2011-11-09 2012-06-20 中兴通讯股份有限公司 移动终端及降低移动终端辐射的方法
CN108027412A (zh) * 2015-07-07 2018-05-11 Q生物公司 场不变定量磁共振特征标志
CN108136195A (zh) * 2015-10-26 2018-06-08 奥海视网膜科技有限公司 使用脉冲的能量源进行生物组织热疗的方法
CN110998971A (zh) * 2017-05-23 2020-04-10 撒加德有限公司 辐射屏蔽罩
WO2022031362A1 (en) * 2020-08-04 2022-02-10 Ojai Retinal Technology, Llc System and method for preventing or treating alzheimer's and other neurodegenerative diseases

Non-Patent Citations (3)

* Cited by examiner, † Cited by third party
Title
吕冬翔，苏东林: "被试品辐射发射试验点位对测试结果影响", 《北京航空航天大学学报》 *
张盛庆宇，李志强，舒逍: "低强度微波辐射对脑组织和脑细胞影响的研究进展", 《职业与健康》 *
王孟均等: "蜂窝移动通信基站电磁辐射对人体影响的探讨", 《中国无线电》 *

Cited By (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN117236110A (zh) * 2023-04-13 2023-12-15 安徽大学 基于电热多物理场随机fdtd的电磁辐射强度评估方法

Also Published As

Publication number Publication date
CN114818412B (zh) 2022-09-27

Similar Documents

Publication Publication Date Title
Kiss et al. Adaptive CAD model (re-) construction with THB-splines
CN114818412B (zh) 一种人脑受电磁辐射的高精度计算方法
Carpio et al. Hybrid topological derivative and gradient-based methods for electrical impedance tomography
Li et al. Influence of the hand on the specific absorption rate in the head
Bruns Three-dimensional simulation and experimental verification of a reverberation chamber
Harmon et al. Adjoint-based accelerated adaptive refinement in frequency domain 3-D finite element method scattering problems
CN112883625A (zh) 一种系统结合的adi-fdtd仿真方法、装置及相关组件
Wu et al. A model-based regularized inverse method for ultrasonic B-scan image reconstruction
Ala et al. A marching-on in time meshless kernel based solver for full-wave electromagnetic simulation
CN110377974A (zh) 雷达波吸收复合板件结构优化方法
Harmon et al. Accelerated adaptive error control and refinement for SIE scattering problems
JP2004038774A (ja) 電磁環境解析方法、電磁環境解析プログラム及びそのプログラムを記録した記録媒体
CN112232001B (zh) 一种集成电路超宽频谐振响应的自适应确定方法及系统
CN105447214A (zh) 器件参数的确定方法和装置
CN114626268B (zh) 一种高精度强电磁脉冲传播过程时域计算方法
Liu et al. Solution of volume–surface integral equation accelerated by MLFMA and skeletonization
Christ et al. Reflection properties of the human skin from 40 to 110 GHz: a confirmation study
CN112347679B (zh) 一种带有电磁色散的超材料仿真方法及装置
WO2011154390A1 (en) Device, method and software progarm for calculating an interaction integral in a conductive medium
CN110162896B (zh) 一种复频移卷积实现完全匹配层的方法
Zhao et al. A hybrid NFM/MoM full-wave analysis of layered prolate head model exposed to handset antenna
CN112147685A (zh) 一种基于波动方程的正演模拟方法及装置
Ali et al. Body area networks at radio frequencies: Creeping waves and antenna analysis
CN113255190B (zh) 一种不规则非均匀介质体雷达散射截面的计算方法
CN116484724B (zh) 基于高阶矩量法电磁建模参数的天线检测方法

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant