突发公共卫生事件下的城市社会经济活动时空变化分析

沈勇杰; 杨春成; 尚海滨; 徐立; 王泽凡

doi:10.13203/j.whugis20220016

突发公共卫生事件下的城市社会经济活动时空变化分析

沈勇杰^1,,
杨春成^{1, 2, 3, ,},
尚海滨¹,
徐立¹,
王泽凡²

1.
中国地质大学（武汉）国家地理信息系统工程技术研究中心,湖北武汉,430078
2.
中国地质大学（武汉）地质探测与评估教育部重点实验室,湖北武汉,430074
3.
中国地质大学（武汉）计算机学院,湖北武汉,430078

基金项目:

国家自然科学基金 42171438

详细信息

作者简介:
沈勇杰,硕士,主要从事夜光遥感与时空数据挖掘等方面的研究。syjjj1016@163.com

通讯作者:
杨春成,博士,研究员。yangcc@cug.edu.cn

计量
- 文章访问数: 10
- HTML全文浏览量: 1
- PDF下载量: 5
出版历程
- 收稿日期: 2023-10-05
- 刊出日期: 2025-03-04

Temporal and Spatial Variation Analysis of Urban Socioeconomic Activity in Public Health Emergencies

1.
National Engineering Research Center of Geographic Information System, China University of Geosciences (Wuhan),Wuhan430078, China
2.
Key Laboratory of Geological Survey and Evaluation of Ministry of Education, China University of Geosciences (Wuhan), Wuhan430074, China
3.
School of Computer Science, China University of Geosciences (Wuhan),Wuhan430078, China

摘要

摘要:
以美国洛杉矶、达拉斯和纽约3个主要城市群为研究区域,基于夜光数据构建平均夜光强度差异值和夜光稳定度两个指标,反映突发公共事件持续期间不同激增期内城市夜光水平的时空变化特征。采用客观赋权法融合移动大数据构建出行意愿指数和场景流动指数,表征区域流动性在不同激增期内的时空差异特征。研究结果表明,城市夜光水平在第一次激增期内降幅较大,城市经济活动受影响更显著,在第二次激增期中,大部分城市经济活动已经恢复至良好水平；城市社会活动变化在第一次激增期内主要表现为场景流动指数的降低,在第二次激增期中,不论是场景流动指数还是出行意愿指数都出现了不同程度的降低。根据激增时间点结合不同地区的社交疏远政策,可以对流动性的变化趋势做出良好的解释。通过夜光数据与移动大数据来表示城市经济社会活动变化,可以了解突发公共卫生事件的反弹对城市社会经济活动的影响,以期为未来可能出现的新型突发公共卫生事件制定对应的社会经济政策提供科学参考。
- 突发公共卫生事件 /
- 夜光稳定度 /
- 场景流动指数 /
- 出行意愿指数 /
- 社会经济活动
Abstract:
Objectives
This paper aims to explore the characteristics of changes in urban socioeconomic activities during a public health emergency. Three major urban agglomerations, including Los Angeles, Dallas and New York, are selected as the study areas.
Methods
First, the difference value of average nighttime light intensity and nighttime light stability are constructed based on nighttime light data to reflect the spatiotemporal variation characteristics of urban nighttime light levels during different surge periods of the duration of public emergencies. Then, based on the objective weighting method and mobility big data, the travel willingness index and scene mobility index are constructed to characterize the spatial and temporal differences of regional mobility in different surge periods.
Results
The results show that, in the first surge period, the level of urban nighttime light decreases greatly and urban economic activities are affected significantly, while in the second surge period, economic activities in most cities recover to a good level. The change of urban social activities in the first surge period is mainly reflected as the decrease of scene mobility index, while in the second surge period, both scene mobility index and travel willingness index decrease to different degrees. Combined with social distancing policies in different regions, the trend of mobility can be well explained according to the surge.
Conclusions
By using nighttime light data and mobility big data to represent the changes of urban social and economic activities, we can understand the impacts of the rebound of public health emergency on urban social and economic activities, in order to provide scientific reference for the development of social and economic policies corresponding to the possible public health emergencies in the future.
- public health emergency /
- nighttime light stability /
- scene mobility index /
- travel willingness index /
- socioeconomic activity
http://ch.whu.edu.cn/cn/article/doi/10.13203/j.whugis20220016

HTML全文

随着海洋强国战略的不断深入,水下潜航器在军事与民事的广泛应用已成为国家海洋发展的战略重点,如何实现水下潜航器长期、自主、隐蔽、高精度航行是满足各项任务要求的先决条件。水下潜航器常使用惯性导航技术进行无源定位导航,但惯性导航系统在长时间导航过程中由于误差累积,定位导航精度大打折扣,重力辅助惯性导航技术具有弥补惯导缺陷、有效修正惯导累积误差的优势,将成为未来水下潜航器定位导航的重点研究方向之一^[1-3]。重力辅助惯性导航的基本原理是利用水下潜航器搭载的相应重力测量仪器实测数据和基准图上惯导指示位置图上数据,通过匹配算法对惯性导航位置和元器件误差进行修正。常用的重力量有重力异常与扰动重力梯度,重力异常匹配算法由于只有一个匹配量,易受到各类误差和载体运动异常的影响。而重力梯度作为重力位的二阶导数,具备6个分量（其中5个独立分量）可以组合搭配进行匹配导航,且重力梯度不易受到水下潜器载体自身加速度的影响,因此扰动重力梯度匹配惯性导航理论上精度更高^[4-5]。

海洋扰动重力梯度基准图是扰动重力梯度辅助惯性导航得以实现的先决工作,如何快速获取高精度扰动重力梯度基准图成为研究重点。本文依据边值问题理论,基于移去-恢复技术将实际扰动重力梯度视为由两部分组成,即地球重力场模型计算得到的模型扰动重力梯度 $T_{α β}^{M}$ 和剩余扰动重力梯度 $δ T_{α β}$ （ $α, β = x, y, z$ ）；然后根据 $T_{α β}^{M}$ 和 $δ T_{α β}$ 的计算特性,使用基于中央处理器（central processing unit,CPU）的OpenMP（open multi-processing）和基于图形处理器（graph processing unit,GPU）的统一计算设备架构（compute unified device architecture,CUDA）并行算法,设计CPU+GPU混合编程并行计算方案,提出数组收缩膨胀方法解决了缔合勒让德函数递推过程内存读写冲突问题,并通过引入Hilbert空间填充曲线将二维格网的基准图降维至一维格网,减少GPU计算负荷和显存读取时间,最终实现全球海洋扰动重力梯度基准图的快速构建。

1 基于CPU的模型扰动重力梯度并行方案

随着计算机技术的快速发展,计算机硬件更新迭代频繁,CPU作为整个计算机系统的计算和控制核心发展极为迅速,已由最初的单核发展到了16核、32核乃至更多核心的众核CPU^[6]。其中OpenMP就是为多核/共享内存处理器设计开发的,是一个基于共享存储器的并行环境^[7]。通常运算程序仅使用CPU的一个线程（即串行程序）,对于硬件快速发展的当下,造成了计算资源的极大浪费。OpenMP则可以通过对串行程序进行简单的改进,实现多线程并行计算任务。同时共享存储系统也有着其天然的缺点,在进行并行计算时,由于存储系统由各个线程共享,因此在计算过程中易导致内存的读写冲突,因此如何对共享内存进行分配和隔离成为需要研究的重点。

1.1 模型扰动重力梯度计算

考虑到梯度的对称性以及主对角线元素之和为0,构建海洋扰动重力梯度基准图通常选取计算6个分量（5个独立分量,1个计算检核条件）,各个分量计算过程与方法相似,本文仅以其中1个分量 $T_{y y}$ 为例介绍扰动重力梯度计算过程。模型扰动重力梯度计算过程可以分为两步：缔合勒让德函数的计算和球谐级数的计算。

模型扰动重力梯度的球谐级数计算形式是较为简单的级数求和形式,可以直接进行并行化处理,其表达式为：

\{\begin{array}{l} T_{y y}^{M} = \sum_{m = 0}^{N} (E_{m}^{y y} c o s (m λ) + F_{m}^{y y} s i n (m λ)) \\ E_{m}^{y y} = \frac{G M}{r^{3}} {\sum_{n = 2}^{N} (\frac{R}{r})}^{n} {\bar{C}}_{n m}^{*} P_{n m}^{y y} (s i n φ) \\ F_{m}^{y y} = \frac{G M}{r^{3}} {\sum_{n = 2}^{N} (\frac{R}{r})}^{n} {\bar{S}}_{n m} P_{n m}^{y y} (s i n φ) \end{array}

(1)

式中,n、m分别表示截断阶数； $R$ 表示平均地球半径； $r$ 表示地心向径,海洋处近似认为 $R / r \approx 1$ ； $φ$ 表示地心纬度； $λ$ 表示地心经度； ${\bar{C}}_{n m}^{*}$ 表示剔除正常重力场偶阶带谐项的模型系数； ${\bar{S}}_{n m}$ 为地球重力场模型系数； $P_{n m}^{y y} (s i n φ)$ 表示计算 $T_{y y}^{M}$ 所使用的勒让德函数。

$P_{n m}^{y y} (s i n φ)$ 的计算需要进行递推运算,即后一步的计算需要前一步的结果作为起算值,因此若直接对其计算过程进行并行化处理,会导致内存读写冲突造成计算崩溃,且其在极点附近（即 $φ$ 接近90°）会产生奇异,因此需要首先对其去奇异,然后解决递推计算并行化问题。 $P_{n m}^{y y} (s i n φ)$ 去奇异后的递推计算式为：

P_{n m}^{y y} (s i n φ) = \{\begin{array}{l} - \frac{1}{2} \sqrt[]{\frac{n (n - 1) (n + 1) (n + 2)}{2}} s i n^{2} φ {\bar{P}}_{n, 2} (s i n φ) - \sqrt[]{\frac{n (n + 1)}{2}} c o s φ s i n φ {\bar{P}}_{n, 1} (s i n φ) - \\ \frac{1}{2} (n + 1) (n s i n^{2} φ + 2) {\bar{P}}_{n, 0} (s i n φ), m = 0 \\ - \frac{1}{4} \sqrt[]{(n - 2) (n - 1) (n + 2) (n + 3)} s i n^{2} φ {\bar{P}}_{n, 3} (s i n φ) - \sqrt[]{(n - 1) (n + 2)} c o s φ s i n φ \cdot \\ {\bar{P}}_{n, 2} (s i n φ) - \frac{1}{4} [(n - 1) s i n^{2} φ + 4] (n + 2) {\bar{P}}_{n, 1} (s i n φ), m = 1 \\ - \frac{1}{4} \sqrt[]{(n - 3) (n - 2) (n + 3) (n + 4)} {\bar{P}}_{n, 4} (s i n φ) - \frac{1}{2} (n^{2} + 3 n + 6) {\bar{P}}_{n, 2} (s i n φ) - \\ \frac{1}{4} \sqrt[]{2 n (n - 1) (n + 1) (n + 2)} {\bar{P}}_{n, 0} (s i n φ), m = 2 \\ - \frac{1}{4} \sqrt[]{(n - m - 1) (n - m) (n + m + 1) (n + m + 2)} {\bar{P}}_{n, m + 2} (s i n φ) - \\ \frac{1}{4} \sqrt[]{(n + m - 1) (n + m) (n - m + 1) (n - m + 2)} {\bar{P}}_{n, m - 2} (s i n φ) - \\ \frac{1}{2} (n^{2} + 3 n + m^{2} + 2) {\bar{P}}_{n, m} (s i n φ), m > 2 \end{array}

(2)

式中, ${\bar{P}}_{n, m} (s i n φ)$ 表示纬度 $φ$ 处的缔合勒让德函数。

由式（2）即可计算全球任意点位的缔合勒让德函数,随后结合CPU单线程计算能力和整体调控能力强的特点,设计利用OpenMP算法对区域模型扰动重力梯度进行并行计算,对缔合勒让德函数递推和球谐级数计算进行整体并行,充分利用CPU计算性能和任务调控能力。

1.2 基于OpenMP的CPU并行方案

由于缔合勒让德函数的计算是一个递推的过程,其计算过程难以进行拆分,需要将整个勒让德函数递推过程作为一个整体,将不同位置点缔合勒让德函数计算作为不同的并行单元（即每一个线程执行的最小单元）,因此每一个并行单元中至少需要划分一部分独立存储空间用来存储缔合勒让德函数,并且并行单元的个数也会直接影响所需存储空间的大小,在计算过程中为保证一定精度,如果缔合勒让德函数计算至2 160阶,共2 160×2 160/2个双精度型的数据需要存储,一个双精度型数据需要占用8 byte的存储空间,此时需要开辟一个17.8 Mbyte的存储空间,这仅仅是一个并行单元只存储一个 ${\bar{P}}_{n m} (s i n φ)$ 所需要的空间^[6]。显然,计算过程不仅仅只需要一个这样的数组,而且并行单元的个数也会因并行方案和计算量的变化而变化,GPU的显存往往没有计算机内存大,且GPU显存通常是不可拓展的^[7],因此勒让德函数计算的并行方案更适合使用CPU并行算法,而后续球谐级数计算部分,为保证计算整体性,减少数据在CPU和GPU之间的传输所消耗的时间,均采用基于OpenMP的CPU并行算法。

确定并行算法后,需要针对计算过程和使用的并行算法进行并行方案设计,方案设计过程中应着重解决内存占用过高和读写冲突问题。为有效解决上述两个问题,本文提出了计算数组收缩膨胀方法如图1所示。

图 1 数组收缩膨胀示意图

Figure 1. Schematic Diagram of the Shrink and Expansion Method of Array

下载: 全尺寸图片幻灯片

以存储缔合勒让德函数 ${\bar{P}}_{n m} (s i n φ)$ 为例,由于缔合勒让德函数所占内存随着截断阶数和级次的增加呈二次方增长,通常计算缔合勒让德函数时,需要开辟一个二维数组P[n][m]对其进行存储,但由缔合勒让德函数特性 $m \leq n$ 可知,P[n][m]仅需要存储下三角部分,其余部分需补0,如果直接使用二维数组进行存储,在计算时有一半的内存被浪费。且在进行并行计算时,为了保证每一个线程所计算的数组独立,需要根据并行单元个数开辟相应个数的数组进行计算。计算过程中还需要对每个线程所用数组进行准确索引,避免因多个线程同时调用一个内存位置而产生内存读写冲突的问题。对此,本文提出的数组收缩膨胀方法利用CPU完成模型扰动重力梯度的计算,既可有效减少内存空间的浪费,又可解决不同线程所用数组的索引问题。具体步骤如下：

1）数组收缩：利用二维矩阵存储 ${\bar{P}}_{n m} (s i n φ)$ ,当阶数为N时,用以存储P[n][m]是一个N×N的二维数组,由于其仅存储下三角部分,一半的内存空间被占用但并未使用,而一维数组则不存在这种问题,因此本文考虑将二维数组的下三角部分以一维数组的形式进行存储,即P[nm],此时一维数组P[nm]仅需占用（N+1）×N/2的内存空间。如图1所示,将二维数组的下三角部分以行为单位,每行内从左至右依次排列到收缩维度后的一维数组中,易得出二维数组P[n][m]中对应位置的勒让德函数在一维数组P[nm]中的对应元素为P[n×（n+1）/2+m]。

2）数组膨胀：在进行收缩后,每一个并行单元都需要独立的P[nm]进行计算,即有t个并行单元就需要t个（N+1）×N/2大小的一维数组,如果简单地开辟t个一维数组P₁[nm],P₂[nm],…,P_t[nm]用来存储,则在进行并行计算过程中不易对其进行索引,易导致因多线程重复使用相同内存位置而产生内存读写冲突。因此本文提出数组膨胀方法,将收缩过后的一维数组膨胀至二维数组P[t][nm],在二维数组中通过索引行数t将不同内存位置分配给不同的线程,既可以解决索引问题,也能避免内存读写冲突问题。

2 基于GPU的剩余扰动重力梯度并行方案

尽管CPU有着强大的单线程计算能力和较大的计算机内存,对于整体并行有着先天优势^[8]。但是受制于材料与工艺,其计算主频和总核心数难以呈量级式增长^[9]。对于单个并行单元计算并不复杂但计算单元众多的任务,CPU的并行加速效果有限,但却十分适合GPU^[10]。其中具有代表性的就是英伟达公司为其旗下GPU产品开发的异构并行平台CUDA^[11]。目前GPU的硬件发展十分迅速,GPU内线程块和线程（CUDA核心）数已由几百提升至上千个,以Tesla V100显卡为例,其CUDA核心数已达到5 120个^[12-14]。并且GPU相较于并行计算机成本较低,可以普遍被大众接受。但是,GPU也存在缺点,尽管GPU有大量的CUDA核心,但是其每一个核心的计算能力无法与CPU相提并论,因此在利用GPU并行计算时,需要对并行单元进行细化,尽可能在充分利用GPU计算能力的条件下,使得CUDA核心运算量不超负荷。在剩余扰动重力梯度的计算过程中,每一个格网点的计算互不相关,可以独立进行计算且单个格网点剩余扰动重力梯度的计算量不大,但格网点数众多,因此不同于CPU的整体并行,GPU则选择计算单一格网的剩余扰动重力梯度作为一个并行单元,增加并行单元个数,减少单个并行单元计算量。

2.1 剩余扰动重力梯度计算

依据移去-恢复技术思路,其目的就是移去扰动重力梯度中的长波分量,而剩余部分则利用计算点周围积分半径内流动点的重力异常剩余量 $δ Δ g$ （ $δ Δ g = Δ g - Δ g_{M}$ , $Δ g$ 表示实际重力异常, $Δ g_{M}$ 表示模型重力异常）进行计算。由于中央区计算存在奇异情况,因此单个并行单元计算过程需要分为内区和中央区。

1）内区剩余扰动重力梯度计算

利用剩余重力异常计算剩余扰动重力梯度的公式为：

δ T_{y y} = \frac{R}{4 π} \iint_{σ}^{} [\frac{1}{r^{2}} (\frac{\partial^{2} S (r, ψ)}{\partial ψ^{2}} {(\frac{\partial ψ}{\partial φ})}^{2} + \frac{\partial S (r, ψ)}{\partial ψ} \frac{\partial^{2} ψ}{\partial φ^{2}}) + \frac{1}{r} \frac{\partial S (r, ψ)}{\partial r}] δ Δ g d σ

(3)

\{\begin{array}{l} \frac{\partial S (r, ψ)}{\partial r} = - \frac{t^{2}}{R} [\frac{1 - t^{2}}{D^{3}} + \frac{4}{D} + 1 - 6 D - t c o s ψ (13 + 6 l n \frac{1 - t c o s ψ + D}{2})] \\ \frac{\partial S (r, ψ)}{\partial ψ} = - t^{2} s i n ψ [\frac{2}{D^{3}} + \frac{6}{D} - 8 - 3 \frac{1 - t c o s ψ - D}{D s i n ψ} - 3 l n \frac{1 - t c o s ψ + D}{2}] \\ \frac{\partial^{2} S (r, ψ)}{\partial ψ^{2}} = - t^{2} c o s ψ [\frac{2}{D^{3}} + \frac{6}{D} - 8 - 3 \frac{1 - t c o s ψ - D}{D s i n^{2} ψ} - 3 l n \frac{1 - t c o s ψ + D}{2}] + \\ t^{3} s i n^{2} ψ [\frac{6}{D^{5}} + \frac{6}{D^{3}} + 3 \frac{D - 1}{D^{2} s i n^{2} ψ} + 3 \frac{D + 1}{D (1 - t c o s ψ + D)} - 3 \frac{1 - t c o s ψ - D}{D s i n^{2} ψ} (\frac{2 c o s ψ}{t s i n^{2} ψ} + \frac{1}{D^{2}})] \end{array}

(4)

式中, $σ$ 表示整个积分区域； $d σ$ 为球面上的流动面元； $ψ$ 是球外计算点与流动面元的球心角距； $S (r, ψ)$ 表示Stokes核函数,其解析式可表示为：

S (r, ψ) = t [\frac{2}{D} + 1 - 3 D - t c o s ψ (5 + 3 l n \frac{1 - t c o s ψ + D}{2})]

(5)

式中, $ψ = c o s^{- 1} (s i n φ s i n φ^{'} + c o s φ c o s φ^{'} c o s (λ^{'} - λ))$ , $φ^{'}$ 、 $λ^{'}$ 分别为球面流动面元的纬度与经度； $t = R / r$ ； $D = \sqrt[]{1 - 2 t c o s ψ + t^{2}}$ 。

利用式（3）与式（4）即可计算全球任意点内区剩余扰动重力梯度 $δ T_{y y}$ 分量。

2）中央区剩余扰动重力梯度计算

由式（4）可知,Stokes核函数及其导数在计算过程中均用到了 $1 / D$ ,因此在流动点接近计算点时, $S (r, ψ)$ 及其导数会异常增大,导致奇异性,这点在计算径向扰动重力梯度时尤为明显^[15-16]。因此,在计算中央区扰动重力梯度时,不能使用原公式计算。本文选用Taylor级数展开法对中央区奇异问题进行求解。

首先,区域大地水准面可以近似看作平面,则大地水准面上的计算点坐标为 $(x, y, 0)$ ,并将中央区视为以计算点为圆心、邻近近1′~2′距离为半径的圆形区域,则半径 $s_{0} = R ψ_{0}$ , $ψ_{0}$ 为中央区圆心与边缘位置的球心角距。在对中央区进行格网化后,半径 $s_{0}$ 可近似等于 $\sqrt[]{Δ x \times Δ y}$ ^[16]。

将直角坐标转化到以计算点为原点的极坐标系下,则流动点在该极坐标系下的坐标为 $(s, α)$ ,其中 $s = \sqrt[]{(x^{'} {- x)}^{2} + (y^{'} {- y)}^{2}}$ , $(x, y)$ 为计算点坐标,（ $x^{'}$ , $y^{'}$ ）为流动点坐标,故易得 $x = s c o s α$ , $y = s s i n α$ , $d x^{'} d y^{'} = s d s d α$ 。因此,可得中央区剩余扰动重力梯度计算式为：

δ T_{y y} = \frac{1}{2 π} \int_{0}^{2 π} \int_{0}^{s_{0}} \frac{δ Δ g}{s^{2}} (3 s i n^{2} α - 1) d s d α

(6)

对式（6）中 $δ Δ g$ 在计算点处进行Taylor级数展开可得：

\begin{matrix} δ Δ g = δ Δ g_{P} + x g_{1} + y g_{2} + \frac{1}{2} (x^{2} g_{11} + \\ 2 x y g_{12} + y^{2} g_{22}) + R_{δ Δ g} \end{matrix}

(7)

式中, $δ Δ g_{P}$ 表示泰勒级数展开点P的重力异常剩余量； $R_{δ Δ g}$ 表示阶段误差；系数参数 $g_{1} = \frac{\partial Δ g}{\partial x}$ , $g_{2} = \frac{\partial Δ g}{\partial y}$ , $g_{11} = \frac{\partial^{2} Δ g}{\partial x^{2}}$ , $g_{22} = \frac{\partial^{2} Δ g}{\partial y^{2}}$ , $g_{12} = \frac{\partial^{2} Δ g}{\partial x \partial y}$ 。

利用上述直角坐标与极坐标转换关系,可将式（7）转换为极坐标形式：

δ Δ g = δ Δ g_{P} + s (g_{1} c o s α + g_{2} s i n α) + \frac{s^{2}}{2} \cdot (g_{11} c o s^{2} α + 2 g_{12} s i n α c o s α + g_{22} s i n^{2} α) + R_{δ Δ g}

(8)

为求得式（8）中系数参数,仅考虑线性项（假设高阶项所占比重较低）,利用剩余重力异常格网数据可得：

\{\begin{array}{l} g_{1_{p, q}} = \frac{δ Δ g_{p + 1, q} - δ Δ g_{p - 1, q}}{2 Δ x} \\ g_{2_{p, q}} = \frac{δ Δ g_{p, q + 1} - δ Δ g_{p, q}}{2 Δ y} \\ g_{11_{p, q}} = \frac{g_{1_{_{p + 1, q}}} - g_{1_{p - 1, q}}}{2 Δ x} \\ g_{12_{p, q}} = \frac{g_{1_{p, q + 1}} - g_{1_{p, q - 1}}}{2 Δ y} \\ g_{22_{p, q}} = \frac{g_{2_{p, q + 1}} - g_{2_{p, q - 1}}}{2 Δ y} \end{array}

(9)

式中,p表示格网行数, $p = 1,2,$ …,P；q表示格网列数, $q = 1,2,$ …,Q。

将式（8）、式（9）代入式（6）可得：

δ T_{y y} = \frac{5 g_{22} - g_{11}}{16} s_{0}

(10)

利用式（10）即可计算全球任意点中央区剩余扰动重力梯度 $δ T_{y y}$ 分量。

2.2 基于CUDA的GPU并行算法

在确定了并行计算单元后,设计利用CUDA算法对剩余扰动重力梯度进行并行计算,在保证GPU单线程计算不超负荷的情况下,充分发挥GPU每一个核心的计算能力。为了区分于CPU程序任务,CUDA有一套独立的核函数形式：kernel_function<<<num_block,num_thread>>>

（param1,param2）,其中 $n u m_b l o c k s$ 表示并行计算时分配的线程块个数； $n u m_t h r e a d s$ 表示每个线程块中所包含的CUDA核心个数。为了更好地控制和调配这些线程块和CUDA核心,CUDA可以对使用的线程块和线程进行编码从而精确划分计算任务。线程块与线程编码如图2所示。

图 2 CUDA线程索引示意图

Figure 2. Diagram of CUDA Thread Index

下载: 全尺寸图片幻灯片

CUDA为了在线程格网中索引特定的线程,根据线程层级划分4个索引变量,其中（ $T h r e a d_{I d, x}$ , $T h r e a d_{I d, y}$ ）分别为线程在线程块中的索引,（ $B l o c k_{I d, x}$ , $B l o c k_{I d, y}$ ）分别为线程所在线程块在线程格网中的索引,D_block,_x和D_block,_y分别表示线程块在x、y方向的维度^[6]。根据上述索引变量,

可以得到任意线程在线程格网中的x、y方向索引计算式为：

\{\begin{array}{l} I d x = B l o c k_{I d, x} \cdot D_{b l o c k, x} + T h r e a d_{I d, x} \\ I d y = B l o c k_{I d, y} \cdot D_{b l o c k, y} + T h r e a d_{I d, y} \end{array}

(11)

图2以一个包含3×5个线程块的线程格网为例,每个线程块中包含3×3个线程。以图2中红色线程为例,其所在的线程块x、y方向维度均为3,其所在线程块在线程格网中x、y方向的索引均为2,该线程在其所在线程块中x、y方向的索引均为1,因此该线程的Idx=7,Idy=7。

在计算过程中由于格网点数众多,利用二维线程格网进行索引增加了线程调度过程的耗时^[6],因此考虑通过一定方法将二维格网进行降维处理,利用一维线程格网进行计算以最大限度地提高计算效率。又由于在进行降维过程中会使原来相邻的格网点分离,使得整体空间聚簇性变差^[17-23],即原相邻格网在一维格网中编码序号相隔较远,因此在存储器中索引时需要读取更多的数据。尤其是在使用GPU时,由于显存读写速度不及主板内存^[6],因此尽量减小二维格网中临近格网在一维格网中的编码差值,提高空间聚簇型,可以有效提高GPU的计算效率。经国内外学者对比研究表明,在众多空间填充曲线中,Hilbert空间填充曲线的空间聚簇性最好^[24-26],二维格网临近点通过Hilbert曲线进行降维后在一维格网中的编码能够尽可能保持接近,有效提高二维格网在一维存储器（显存）上的索引效率^[27]。因此,本文考虑使用Hilbert空间填充曲线对二维格网进行降维处理,其空间填充方式如图3所示。

图 3 Hilbert空间填充曲线示意图

Figure 3. Schematic Diagram of Hilbert Space Filling Curve

下载: 全尺寸图片幻灯片

Hilbert空间填充曲线建立二维以及高维空间格网和一维格网间的映射关系,这样的映射关系序列码被称为Hilbert空间排列码^[23-25]。常用的二维Hilbert空间填充曲线有三阶,如图3所示。1阶Hilbert曲线由4个相邻格网构成,空间曲线不重复地贯穿整个格网；2阶Hilbert曲线由4个1阶曲线通过旋转构成,曲线贯穿区域；同样,3阶Hilbert曲线是由1个2阶曲线通过旋转构成,曲线贯穿整个区域。为了对计算区域进行有效降维编码,避免多余的计算量,本文仅使用最高为3阶的Hilbert曲线,通过将计算区域划分为K个3阶Hilbert曲线首尾相连的方式对计算区域进行降维。Hilbert空间填充曲线的格网编码与解码可采用递归算法和迭代算法,因此,只需要对3阶Hilbert曲线格网进行顺序编码,得到待求格网点所在的3阶Hilbert格网编码以及其在该格网中的编码,即可在一维数组中对其进行索引。以图3中红色格网点为例,其在3阶Hilbert曲线上的值为53。经过上述并行方案,即可利用GPU完成剩余扰动重力梯度的计算。

3 实验验证与分析

为验证本文所提方法的有效性,在下述硬件、软件实验平台对该算法进行实验验证：实验环境为高性能台式机,处理器为Intel Xeon（R） Gold 6130 @ 2.10 GHz；GPU为Tesla V100,内存64 GB；软件环境为64位Windows 10操作系统,VS 2013+CUDA 9.2平台。

3.1 精度分析

为验证本文所提并行算法精度,选取6°×6°范围研究区（海洋格网范围为8°×8°,位于6°N~14°N,158°E~166°E）。该区域海深最深6 872 m,最浅0 m,区域平均高程-3 657.22 m；该区域重力异常最大值313.8 mGal,最小值-59.3 mGal,平均重力异常12.4 mGal。区域重力异常与海深如图4所示。

图 4 实验区域物理特性图

Figure 4. Physical Characteristic Map of the Study Area

下载: 全尺寸图片幻灯片

目前,GPU在进行单精度浮点运算的效率要远远高于双精度浮点运算^[6],对于同一区域扰动重力梯度基准图的计算,使用单精度数据计算剩余扰动梯度的耗时仅为利用双精度数据计算耗时的1/4,但由于在计算剩余扰动重力梯度过程中存在大量正、余弦、正切、开方等函数,仅使用单精度浮点计算结果不能满足计算精度的要求,有限的有效位数带来的舍入误差会变得更加明显,甚至影响计算精度。为此,本文利用GPU分别使用单精度浮点数据和双精度数据,计算实验区域剩余扰动重力梯度6个分量,并利用CPU使用双精度数据串行计算实验区域剩余扰动重力梯度,以此结果作为真值,与GPU计算结果作差并进行对比分析,精度分析结果见表1（1 E=10^-9 s^-2）。

表 1 不同数据类型并行计算精度分析/E

Table 1. Accuracy Analysis of Parallel Computing for Different Data Types/E

数据类型	最大差值	最小差值	平均差值	RMS
双精度数据	1×10^-6	-1×10^-6	<1×10^-12	<1×10^-12
单精度数据	5.32	-1.36	0.18	0.42

下载: 导出CSV

| 显示表格

为进一步验证并行方案数学模型的可靠性和正确性,根据上述计算结果,基于 $T_{x x} + T_{y y} + T_{z z} = 0$ 的理论约束条件,分别将单精度浮点数据与双精度浮点数据运算结果的内符合情况进行统计,结果如图5所示,数值特征统计见表2。

图 5 不同数据类型计算

T_{x x} + T_{y y} + T_{z z}

结果

Figure 5. Result of

T_{x x} + T_{y y} + T_{z z}

for Different Data Types

下载: 全尺寸图片幻灯片

表 2

T_{x x} + T_{y y} + T_{z z}

数值特征分析/E

Table 2. Analysis of

T_{x x} + T_{y y} + T_{z z}

Numerical Characteristics/E

数据类型	最大值	最小值	平均值	标准差	RMS
单精度数据	4.43	-3.21	0.20	0.023	0.047
双精度数据	0.047	-0.094	<1×10^-4	0.002 6	0.003 9

下载: 导出CSV

| 显示表格

由图5（b）、图5（c）和表2可知,利用双精度数据计算所得 $T_{x x} + T_{y y} + T_{z z}$ 绝对值最大值小于0.1 E,平均值小于1×10^-4 E,标准差为0.002 6 E,均方根值为0.003 9 E, $T_{x x} + T_{y y} + T_{z z}$ 围绕着0随机波动,符合扰动重力梯度的内符合条件,能够在一定程度上表明上述计算方法数学模型的可靠性和软件编程的正确性。由表1、表2和图5（a）可知,利用单精度数据进行计算,由于数据有效位数的限制,由舍入误差带来的剩余扰动重力梯度的误差绝对值最大可达5.32 E,平均误差也有0.18 E,且在利用 $T_{x x} + T_{y y} + T_{z z}$ 进行内符合精度分析时,最小达到-3.21 E,最大达到4.43 E,无法满足海洋扰动重力梯度基准图的精度需求^[4]。而使用双精度数据进行计算,并行计算绝对误差绝对值最大仅为1×10^-6 E,平均误差小于1×10^-12 E,并行计算精度满足海洋扰动重力梯度基准图的计算精度需求。因此在利用GPU进行剩余扰动梯度计算时,应牺牲一定计算效率以满足精度需求,选取双精度数据进行计算。

3.2 效率分析

3.2.1 不同并行方案效率分析

在使用双精度数据保证计算精度的前提下,为验证本文所提并行算法的效率,分别使用串行计算方案、纯CPU并行方案、纯GPU并行方案及CPU+GPU混合并行方案,对实验区域扰动重力梯度进行计算,模型扰动重力梯度计算采用EIGEN6C4地球重力场模型,截断阶数计算至2 160阶^[28]；使用卫星测高数据反演所得DTU18模型数据（1′分辨率）作为实际重力异常值进行剩余扰动重力梯度计算。计算剩余扰动重力梯度时,内区选取积分半径1°范围进行计算,中央区选取中心格网点和围绕其周围的8个格网点共计9个格网区域。进行3次计算,计算耗时取3次结果的平均值。计算结果如图6所示,计算耗时统计见表3。

图 6 实验区域重力梯度图

Figure 6. Disturbing Gravity Gradient Map of the Study Area

下载: 全尺寸图片幻灯片

表 3 不同并行方案计算耗时

Table 3. Time Consuming of Different Parallel Schemes

方案	模型扰动重力梯度计算耗时/s	模型扰动重力梯度计算加速比	剩余扰动重力梯度计算耗时/s	剩余扰动重力梯度计算加速比	总耗时/s	总加速比
串行计算	544.99	1.00	1 582.560	1.00	2 127.55	1.00
纯CPU并行	14.36	37.95	41.270	38.34	55.63	38.24
纯GPU并行	31.63	17.23	1.484	1 066.21	33.11	64.26
CPU+GPU混合并行	14.36	37.95	1.484	1 066.21	15.84	134.26

下载: 导出CSV

| 显示表格

图6显示,利用双精度数据进行计算,4种计算方案的计算结果与该区域重力异常图和海底地形图有相关性,计算精度一致,任意两种计算方案间的最大偏差仅为1×10^-6 E,且均符合 $T_{x x} + T_{y y} + T_{z z} = 0$ 的内符合条件,因此4种方案计算精度均可靠。由表3可知,利用传统串行方法计算该区域扰动重力梯度全张量需要2 127.55 s（35.46 min）,利用并行方案进行计算耗时均小于1 min,其中,纯CPU并行方案计算耗时55.63 s,总加速比约38倍；纯GPU并行方案计算耗时33.11 s,总加速比约64倍；CPU+GPU混合并行方案计算耗时15.84 s,总加速比约134倍,加速效果优于单系统并行方案。在扰动重力梯度计算过程中,模型扰动重力梯度部分计算耗时占比25.6%,剩余扰动重力梯度部分计算耗时占比74.4%,利用CPU并行方案计算模型部分和剩余部分的加速比分别为37.95倍和38.34倍,而利用GPU计算两部分的加速比分别为17.23倍和1 066.21倍,因此,本文所提CPU+GPU混合并行方案可以充分利用计算设备的运算能力,取长补短,最大限度地提高扰动重力梯度基准图的构建效率。

3.2.2 不同分辨率计算效率分析

为进一步验证本文算法的有效性,使用不同分辨率数据进行数值实验,计算范围选取2°×2°海洋区域,模型扰动重力梯度计算至2 160阶。计算剩余扰动重力梯度时,内区积分半径为1°,中央区选取中心格网点和围绕其周围的8个格网点共计9个格网区域,分辨率分别为1′、30″、10″,计算耗时如表4所示。

表 4 不同分辨率计算耗时

Table 4. Calculation Time of Different Resolutions

分辨率	串行计算耗时	并行计算耗时	加速比
1′	4.02 min	1.92 s	125.36
30″	52.16 min	22.65 s	138.16
10″	64.85 h	26.67 min	145.86

下载: 导出CSV

| 显示表格

由表4可以看出,随着分辨率的提高,计算耗时呈指数上涨,大范围高分辨率海洋扰动重力梯度基准图难以在短时间内得到,以2°×2°范围海洋扰动重力梯度基准图的构建任务为例进行实验。利用本文所提CPU+GPU并行算法计算不同分辨率海洋扰动重力梯度基准图的加速比均可达到125倍以上,最高达到145倍。完成1′分辨率海洋扰动重力梯度基准图构建仅需1.92 s；完成30″分辨率基准图构建仅需22.65 s,而传统串行方法则需要约52 min；完成10″分辨率基准图构建需要26.67 min,而串行方法则需要2 d以上。证明了本文所提并行方案对不同分辨率扰动重力梯度基准图的构建有着普遍有效性,可以有效提高各个分辨率扰动重力梯度基准图的构建速度。

4 全球海洋扰动重力梯度基准图构建与结果分析

通过上述CPU+GPU混合并行方案,可以实现大范围海洋扰动重力梯度基准图快速构建任务,本文对全球范围、1′分辨率海洋区域的扰动重力梯度基准图进行构建,计算使用EIGEN6C4地球重力场模型,截断阶数计算至2 160阶；使用卫星测高数据反演所得DTU18模型数据（1′分辨率）作为实际重力异常值进行剩余扰动重力梯度计算。计算剩余扰动重力梯度时,内区选取积分半径1°范围进行计算,中央区选取中心格网点和围绕其周围的8个格网点共计9个格网区域。并行算法耗时4.52 h,串行计算耗时超过28 d。计算结果如表5、图7所示。

表 5 全球海洋扰动重力梯度数值特征

Table 5. Numerical Characteristics of Global Ocean Gravity Gradient

重力梯度	最大值/E	最小值/E	平均值/E	中位数/E
T_xx	127.33	-116.93	<1×10^-6	0.082 7
T_yy	100.73	-105.30	<1×10^-6	0.060 7
T_zz	228.28	-114.65	<1×10^-6	-0.020 2
T_xy	80.22	-71.67	<1×10^-6	-0.002 0
T_zx	106.48	-112.24	<1×10^-6	0.001 5
T_zy	106.13	-104.44	<1×10^-6	-0.005 9

下载: 导出CSV

| 显示表格

图 7 全球海洋扰动重力梯度图

Figure 7. Global Ocean Disturbing Gravity Gradient Map

下载: 全尺寸图片幻灯片

由表5、图7可知,全球海洋扰动重力梯度平均值均小于1×10^-6 E,且中位数在0 E附近,说明海洋区域大部分点位扰动重力梯度绝对值较小；扰动重力梯度与海底地形呈相关性,地形起伏剧烈区域,扰动重力梯度值变化较为明显。

5 结语

本文利用CPU+GPU混合并行方法,针对扰动重力梯度计算各过程特点,提出了一种高效的并行计算方案,可以快速高效地完成大范围乃至全球范围海洋扰动重力梯度基准图的构建任务。经实验计算分析表明：

1）本文所提方法通过数组收缩膨胀方法和引入Hilbert空间填充曲线,有效克服了缔合勒让德函数递推计算的并行化问题以及GPU显存快速索引问题,可以充分利用设备计算能力。

2）本文所提方法的并行计算误差小于1×10^-6 E,均方根误差小于1×10^-12 E,并行计算过程产生的误差较小,计算精度较高。所得 $T_{x x} +$

$T_{y y} + T_{z z}$ 绝对值最大值小于0.1 E,平均值小于1×10^-6 E,标准差为0.002 6 E,均方根误差为0.003 9 E, $T_{x x} + T_{y y} + T_{z z}$ 围绕着0随机波动,符合扰动重力梯度的内符合条件,表明该并行方案具有一定的可靠性。

3）利用该并行方案计算6°×6°范围、1′分辨率海洋扰动重力梯度基准图仅需15.84 s,而传统串行方法需要35 min以上；计算2°×2°范围、30″分辨率海洋扰动重力梯度基准图仅需22.65 s,而串行方法则需要约52 min；计算2°×2°范围、10″分辨率基准图需要26.67 min,而串行方法则需要2 d以上时间。利用该CPU+GPU混合并行方案进行海洋扰动重力梯度基准图构建,相较于传统串行计算,加速比最高可达140倍以上,有效减少了计算耗时,提高了计算效率。

4）全球全张量海洋扰动重力梯度基准图可在4.6 h内完成构建,相较于串行算法效率提高了150倍以上,由实验结果可知,全球海洋扰动重力梯度平均值均小于1×10^-6E,中位数在0 E附近。

http://ch.whu.edu.cn/cn/article/doi/10.13203/j.whugis20220016

图 1 3个城市群的每周确诊人数变化曲线

Figure 1. Weekly Change Curves of Confirmed Cases in Three Urban Agglomerations

下载: 全尺寸图片幻灯片

图 2 场景访问指数的构建

Figure 2. Construction of Scene Visiting Index

下载: 全尺寸图片幻灯片

图 3 郡级行政区域的平均夜光强度差异值

Figure 3. DANL in County Areas

下载: 全尺寸图片幻灯片

图 4 不同激增期内不同夜光评价指标的变化幅度

Figure 4. Variation Range of Different Evaluation Indexes in Different Surge Periods

下载: 全尺寸图片幻灯片

图 5 3个城市群的夜光变化度空间分布图

Figure 5. Spatial Distribution of Nighttime Light Intensity Variation in Three Urban Agglomerations

下载: 全尺寸图片幻灯片

图 6 研究区域内SMI与TWI曲线

Figure 6. SMI and TWI Curves of the Study Area

下载: 全尺寸图片幻灯片

图 7 不同激增期内的3个城市群的SMI与TWI曲线

Figure 7. SMI and TWI Curves in Different Surge Periods in Three Urban Agglomerations

下载: 全尺寸图片幻灯片

表 1 3个城市群的激增期时间

Table 1 Surge Time in Three Urban Agglomerations

城市群	第一次激增期	第二次激增期
洛杉矶	2020-06-01—2020-09-13	2020-11-02—2021-02-28
达拉斯	2020-06-15—2020-08-23	2020-11-02—2021-02-28
纽约	2020-03-16—2020-05-10	2020-11-09—2021-02-28

下载: 导出CSV

表 2 3个城市群的夜光数据研究时期

Table 2 Study Period of Nighttime Light Data in Three Urban Agglomerations

城市群	第一次激增期		第二次激增期
城市群	激增开始前期	快速激增时期	激增开始前期	快速激增时期
洛杉矶	2020年4月	2020年7月	2020年9月	2020年12月
达拉斯	2020年5月	2020年7月	2020年10月	2021年1月
纽约	2020年2月	2020年3月	2020年9月	2021年1月

下载: 导出CSV

表 3 基于人口比例的3个城市群移动数据的融合权重

Table 3 Population-Based Fusion Weights of Mobility Data in Three Urban Agglomerations

城市群	行政单位	2020年人口数/人	融合权重
洛杉矶	Los Angeles	10 114 218	0.760 7
洛杉矶	Orange	3 181 516	0.239 3
达拉斯	Collin	1 236 538	0.183 1
	Dallas	2 520 155	0.040 7
	Denton	1 007 795	0.128 2
	Tarrant	2 251 621	0.077 3
纽约	Bronx	1 435 322	0.176 2
	Kings	2 542 010	0.359 2
	Nassau	1 334 702	0.143 6
	New York	1 630 542	0.320 9
	Queens	2 227 117	0.118 0
	Richmond	494 553	0.209 0
	Suffolk	1 559 886	0.109 7
	Westchester	939 762	0.134 0

下载: 导出CSV

表 4 不同激增期的平均夜光强度差异值/(nW·cm^-2·sr)

Table 4 DANL of Different Surge Periods/(nW·cm^-2·sr)

城市群	2019年		2020年
城市群	第一次激增期同期	第二次激增期同期	第一次激增期	第二次激增期
洛杉矶	-0.167 6	3.367 0	-7.662 8	6.591 9
达拉斯	1.844 1	-2.097 0	-0.900 6	6.098 1
纽约	3.448 7	4.782 1	-1.149 1	8.984 0

下载: 导出CSV

表 5 各郡县在不同激增期内的夜光稳定度变化/%

Table 5 Changes of Nighttime Light Intensity Stability in Different Surge Periods in Different Counties/%

城市群	包含行政区域	第一次激增期	第二次激增期
洛杉矶	Los Angeles	23.24	56.68
洛杉矶	Orange	21.51	73.62
达拉斯	Collin	25.06	81.95
	Dallas	33.89	85.36
	Denton	37.42	86.23
	Tarrant	46.05	79.03
纽约	Bronx	16.39	84.23
	Kings	36.63	89.46
	Nassau	13.90	69.29
	New York	59.90	83.08
	Queens	19.63	83.53
	Richmond	6.98	85.31
	Suffolk	20.06	64.60
	Westchester	10.39	46.09

下载: 导出CSV

表 6 基于客观赋权法的移动数据融合权重分配

Table 6 Weight Allocation of Mobility Data Fusion Based on Objective Weighting Method

数据库	分类指标	洛杉矶	达拉斯	纽约
Google	零售店与休闲场所	0.319 3	0.256 4	0.303 9
	食品杂货店与药房	0.367 3	0.402 1	0.384 9
	工作场所	0.313 5	0.341 5	0.311 2
Apple	驾驶	0.269 8	0.396 8	0.291 5
	公共交通	0.446 8	0.306 9	0.441 6
	步行	0.283 3	0.296 3	0.266 9

下载: 导出CSV

参考文献(18)

[1]	GAO S, RAO J M, KANG Y H, et al. Mapping County-Level Mobility Pattern Changes in the United States in Response to COVID-19[J]. SIGSPATIAL Special, 2020, 12（1）： 16-26.
[2]	KRAEMER M U G, YANG C H, GUTIERREZ B, et al. The Effect of Human Mobility and Control Measures on the COVID-19 Epidemic in China[J]. medRxiv, 2020： 2020.03.02.20026708.
[3]	MAIER B F, BROCKMANN D. Effective Containment Explains Subexponential Growth in Recent Confirmed COVID-19 Cases in China[J]. Science, 2020, 368（6492）： 742-746.
[4]	HUANG X, LI Z L, JIANG Y Q, et al. Twitter Reveals Human Mobility Dynamics During the COVID-19 Pandemic[J]. PLoS One, 2020, 15（11）： e0241957.
[5]	钟雷洋, 周颖, 高松, 等. 突发公共卫生事件下的人口流动模式变化识别[J]. 武汉大学学报（信息科学版）, 2024, 49（7）： 1237-1249. ZHONG Leiyang, ZHOU Ying, GAO Song, et al. Identifying Human Mobility Patterns Changes During Public Health Emergencies[J]. Geomatics and Information Science of Wuhan University, 2024, 49（7）： 1237-1249.
[6]	李德仁, 邵振峰, 于文博, 等. 基于时空位置大数据的公共疫情防控服务让城市更智慧[J]. 武汉大学学报（信息科学版）, 2020, 45（4）： 475-487. LI Deren, SHAO Zhenfeng, YU Wenbo, et al. Public Epidemic Prevention and Control Services Based on Big Data of Spatiotemporal Location Make Cities more Smart[J]. Geomatics and Information Science of Wuhan University, 2020, 45（4）： 475-487.
[7]	YILMAZKUDAY H. COVID-19 Spread and Inter-County Travel： Daily Evidence from the U.S[J]. Transportation Research Interdisciplinary Perspectives, 2020, 8： 100244.
[8]	ZHU D S, MISHRA S R, HAN X K, et al. Social Distancing in Latin America During the COVID-19 Pandemic： An Analysis Using the Stringency Index and Google Community Mobility Reports[J]. Journal of Travel Medicine, 2020, 27（8）： taaa125.
[9]	李德仁, 李熙. 论夜光遥感数据挖掘[J]. 测绘学报, 2015, 44（6）： 591-601. LI Deren, LI Xi. An Overview on Data Mining of Nighttime Light Remote Sensing[J]. Acta Geodaetica et Cartographica Sinica, 2015, 44（6）： 591-601.
[10]	LAN T, SHAO G F, TANG L N, et al. Quantifying Spatiotemporal Changes in Human Activities Induced by COVID-19 Pandemic Using Daily Nighttime Light Data[J]. IEEE Journal of Selected Topics in Applied Earth Observations and Remote Sensing, 2021, 14： 2740-2753.
[11]	GHOSH T, ELVIDGE C D, HSU F C, et al. The Dimming of Lights in India During the COVID-19 Pandemic[J]. Remote Sensing, 2020, 12（20）： 3289.
[12]	ELVIDGE C, GHOSH T, HSU F C, et al. The Dimming of Lights in China During the COVID-19 Pandemic[J]. Remote Sensing, 2020, 12（17）： 2851.
[13]	YIN R Y, HE G J, JIANG W, et al. Night-Time Light Imagery Reveals China’s City Activity During the COVID-19 Pandemic Period in Early 2020[J]. IEEE Journal of Selected Topics in Applied Earth Observations and Remote Sensing, 2020, 14： 5111-5122.
[14]	BEYER R C M, FRANCO-BEDOYA S, GALDO V. Examining the Economic Impact of COVID-19 in India Through Daily Electricity Consumption and Nighttime Light Intensity[J]. World Development, 2021, 140： 105287.
[15]	SMALL C, SOUSA D. Spatiotemporal Evolution of COVID-19 Infection and Detection Within Night Light Networks： Comparative Analysis of USA and China[J]. Applied Network Science, 2021, 6（1）： 10.
[16]	STRAKA W, KONDRAGUNTA S, WEI Z G, et al. Examining the Economic and Environmental Impacts of COVID-19 Using Earth Observation Data[J]. Remote Sensing, 2021, 13（1）： 5.
[17]	陈斌, 徐尚昭, 周阳阳, 等. POI与NPP/VIIRS夜光数据空间耦合关系下的城市空间结构分析：以武汉市主城区为例[J]. 测绘通报, 2020（7）： 70-75. CHEN Bin, XU Shangzhao, ZHOU Yangyang, et al. Analysis of Urban Spatial Structure Based on Spatial Coupling Between POI and Nighttime Data： A Case Study of Main Urban Region in Wuhan[J]. Bulletin of Surveying and Mapping, 2020（7）： 70-75.
[18]	宋冬梅, 刘春晓, 沈晨, 等. 基于主客观赋权法的多目标多属性决策方法[J]. 山东大学学报（工学版）, 2015, 45（4）： 1-9. SONG Dongmei, LIU Chunxiao, SHEN Chen, et al. Multiple Objective and Attribute Decision Making Based on the Subjective and Objective Weighting[J]. Journal of Shandong University （Engineering Science）, 2015, 45（4）： 1-9.

施引文献

资源附件(0)

图(7) / 表(6)

计量

文章访问数: 10
HTML全文浏览量: 1
PDF下载量: 5
被引次数: 0

1 基于CPU的模型扰动重力梯度并行方案
1.1 模型扰动重力梯度计算
1.2 基于OpenMP的CPU并行方案
2 基于GPU的剩余扰动重力梯度并行方案
2.1 剩余扰动重力梯度计算
2.2 基于CUDA的GPU并行算法
3 实验验证与分析
3.1 精度分析
3.2 效率分析
3.2.1 不同并行方案效率分析
3.2.2 不同分辨率计算效率分析
4 全球海洋扰动重力梯度基准图构建与结果分析
5 结语

1 基于CPU的模型扰动重力梯度并行方案
1.1 模型扰动重力梯度计算
1.2 基于OpenMP的CPU并行方案
2 基于GPU的剩余扰动重力梯度并行方案
2.1 剩余扰动重力梯度计算
2.2 基于CUDA的GPU并行算法
3 实验验证与分析
3.1 精度分析
3.2 效率分析
3.2.1 不同并行方案效率分析
3.2.2 不同分辨率计算效率分析
4 全球海洋扰动重力梯度基准图构建与结果分析
5 结语

参考文献(18)

施引文献

资源附件(0)

突发公共卫生事件下的城市社会经济活动时空变化分析

作者简介: 沈勇杰,硕士,主要从事夜光遥感与时空数据挖掘等方面的研究。syjjj1016@163.com

通讯作者: 杨春成,博士,研究员。yangcc@cug.edu.cn

计量

出版历程

Temporal and Spatial Variation Analysis of Urban Socioeconomic Activity in Public Health Emergencies

1 基于CPU的模型扰动重力梯度并行方案

1.1 模型扰动重力梯度计算

1.2 基于OpenMP的CPU并行方案

2 基于GPU的剩余扰动重力梯度并行方案

2.1 剩余扰动重力梯度计算

2.2 基于CUDA的GPU并行算法

3 实验验证与分析

3.1 精度分析

3.2 效率分析

3.2.1 不同并行方案效率分析

3.2.2 不同分辨率计算效率分析

4 全球海洋扰动重力梯度基准图构建与结果分析

5 结语

计量

出版历程

目录

1 基于CPU的模型扰动重力梯度并行方案

1.1 模型扰动重力梯度计算

1.2 基于OpenMP的CPU并行方案

2 基于GPU的剩余扰动重力梯度并行方案

2.1 剩余扰动重力梯度计算

2.2 基于CUDA的GPU并行算法

3 实验验证与分析

3.1 精度分析

3.2 效率分析

3.2.1 不同并行方案效率分析

3.2.2 不同分辨率计算效率分析

4 全球海洋扰动重力梯度基准图构建与结果分析

5 结语

作者简介:
沈勇杰,硕士,主要从事夜光遥感与时空数据挖掘等方面的研究。syjjj1016@163.com

通讯作者:
杨春成,博士,研究员。yangcc@cug.edu.cn