Laplace

Laplace（共7篇）

Laplace 篇1

1 概述

拉普拉斯变换是工程数学常用的一种积分变换,在工程技术和理论研究等诸多方面有较广泛的应用,特别是对于求解常系数线性偏微分方程的初边值问题。本文基于离散化的思想和蒙特卡罗算法对某一常见函数进行Laplace变换,利用两种并行模式改造并精选出最优的并行计算处理方法。

蒙特卡罗方法是一种随机模拟方法,它以统计理论方法和概率为基础,是使用随机数或伪随机数来解决实际计算问题的方法。Open MP是一种面向共享内存及分布式共享内存的多处理器多线程并行编程语言,是能用于显示指导多线程、共享内存并行的应用程序编程接口,并利用封装好的函数来进行并行计算;Win API并行模式是通过编写线程函数,利用Windows系统中的API接口进行多线程调度来实现并行。通过对数值实验的串并行算法的计算速度,加速效果进行比较,得出并行优化后的计算方法。并进一步考虑将计算中高效的并行方法应用到其他数值问题。

2 拉普拉斯变换描述

定义设f(t)是时间t的函数,当t>0时,对定义在区间(0,+∞)上的实自变量t的函数f(x),乘以e-pt(其中p为复数),然后对t由0到+∞积分,此广义积分若收敛,称函数F(p)

为函数f(t)拉普拉斯变换,记作L[f(t)],即

3 某问题的串并行算法描述

当f(x)=sin(x)时,计算其Laplace变换:

算法改进:

由于e-px的变化速度太快,我们没有办法直接求sin(x)的Laplace变换,故将其离散化形成级数形式,然后利用梯形公式法则的数值积分方法来估算它的值。为了提高精度把积分区间[0,200000]分成若干等份,分点为xn+1=xn+step,形成若干个自区间,再在每个子区间上采用梯形公式。

串行程序:

3.1 API实现

首先在进程内创建线程,它是CPU调度和分配的基本单位。利用WINAPI定义线程函数,在线程函数内分配该线程所要进行的任务,然后将线程函数导入到创建好的线程中运行,计算机通过API接口,针对创建的线程数目,调度相同数目的CPU进行计算,最后将各线程的计算结果合并得到最终结果。

程序中使用多线程时,很多情况下是一些线程进行某些处理操作,而其他线程必须对其处理结果进行了解,若不采取适当措施,其他线程会在线程处理任务结束前就去访问处理结果,因此产生了数据竞争。为避免数据竞争,我们使用临界区进行数据同步,临界区独占对某些共享资源的访问,它只允许一个线程对共享资源进行访问。若有多个线程想同时访问临界区,那么在一个线程进入临界区后其他线程将被挂起。在本问题的实现过程中,将求和运算拖入临界区,来避免频繁介入临界区。

API关键代码如下:

3.2 Open MP实现

利用parallel for循环并行化

采用工作分配的执行方式,利用parallel for循环语句将总的工作量按照一定的方式分配到各个线程,最后将所有线程的结果汇总。Parellel for是对一个完整的for循环进行分割,它根据线程数大小按工作量平均分配,每个线程的工作量为总的工作量/线程数。为了避免产生数据竞争,采用的方法是对竞争变量私有化,添加私有量private(n),再对执行完的程序利用reducation语句进行归约操作,用来更好的收集任务结果。

关键代码

4 蒙特卡罗算法求解

4.1 蒙特卡罗算法

蒙特卡罗法计算多重积分时,采用均匀随机数法得到随机变量是常用方法之一。选取f(p)的方法是取Vs上的均匀分布,即

Vs表示积分区域的体积,此时g(p)=VsG(p)。采用均匀分布随机数计算多重积分时,有。

设区域D:a≤x≤b,c≤y≤d其面积A=(b-a)(d-c),f(x,y)设为区域上的有界函数,不妨设(xi,yi),i=1,2...N为落在区域D上的均匀分布随机数列,用均匀随机数计算二重积分时,当有N充分大时,有。将上式推广到多重积分,有

4.2 基于Laplace变化的蒙特卡洛串行算法描述

Laplace变化格式为单重积分,所以我们在求解Laplace变换过程中我们可以用蒙特卡洛求解多重积分的方法对它求解。根据上述蒙特卡洛算法将Laplace变化最终转化为。然后通过编程求解。

关键串行代码如下:

4.3 基于Open Mp并行算法实现

我们在串行程序的基础上通过Open MP进行并行处理。求解该问题需要计算N个值,然后对N个值进行求和。在计算N个值得过程中数据之间没有依赖关系,所以我们采用四个线程将任务进行分配。在计算每个值之前我们都需要随机产生一个介于积分区间的随机数,为了减少共享内存带来的时间消耗,我们将这个随机数设定为每个线程的私有变量private(x)。然后利用规约reduction(+:s)将N个值进行求和运算,得到最终结果。

关键代码:

5 分析与结论

四线程并行计算:

上述结果可以看出利用离散化方法所花费的时间整体比较少,在串行的基础上对它进行并行实现,使得效率进一步提高。Open MP的加速比达到2.5396,API的加速比达到3.113。利用蒙特卡洛算法计算过程中所用时间比较长。主要因为蒙特卡洛算法主要针对的是多重积分的计算,而该问题是一重积分的运算所以效率较低。但蒙特卡洛算法的串行和并行算法相比较可以看出并行算法效率较高,加速比达到2.498。综合上述结果可以看出基于离散方法的API并行效率最高。

6 结语

对于Laplace变化,可以看出基于离散算法的API的并行效率最高,蒙特卡洛Open MP效率最低。Open MP编程模型比较适合迭代的并行计算。类似的我们还可以通过并行手段解决其他的并行计算问题。

摘要：Laplace变换是一种解决常系数偏微分方程初边值问题的常用方法,在此,利用离散化方法和蒙特卡罗算法求解Laplace变换,运用多线程技术,结合高性能并行计算中Win API,Open MP的并行模式获得最优计算方法。通过数值结果可以看出,蒙特卡罗Open MP方法速度最慢,效率最低,基于离散算的Win API方法计算速度最快,并行效率最高,最终获得结论的并行计算模式应用到其他初值问题。

关键词：Laplace变换,OpenMP,WinAPI,蒙特卡罗算法

参考文献

[1]朱建伟,刘荣.多线程并行快速求解e值得六种方法[J].现代计算机,2013(6).

[2]薛小超,王克,朱朋海.基于多核并行的非线性方程蒙特卡洛计算方法[J].电脑知识与技术,2014(7).

[3]刘辉玲,叶峰.计算多重积分的均匀随机数蒙特卡罗法的实现[J]电脑知识与技术,2008(12).

[4]张继伟.基于MPI的并行拉普拉斯变换算法及其应用研究[D].湖北大学,2011.

Laplace 篇2

Henstock-Kurzweil可积函数的Laplace变换的反演定理

该文建立了Henstock-Kurzweil可积函数的`Laplace变换,讨论了其基本性质及解析性质,得到Henstock-Kurzweil可积意义下的反演公式,并给出反例说明这一结果不能改进.

作 者：苏峰 彭志刚 Su Feng Peng Zhigang 作者单位：湖北大学数学与计算机科学学院,武汉,430062刊 名：数学物理学报 ISTIC PKU英文刊名：ACTA MATHEMATICA SCIENTIA年，卷(期)：27(6)分类号：O177.6关键词：Henstock-Kurzweil可积函数 Laplace变换 反演变换.

循环图的Laplace的谱半径 篇3

设G是一个n阶连通简单图.G的邻接矩阵记为A (G) , 它是一个n×n矩阵, 记作

$A(G)=(a{}_{ij})n\times n‚a{}_{ij}=\{\begin{array}{l}1‚\\ 0.\end{array}A(G)$

的特征值又称为C的特征值;由于A (G) 是非负不可约对称矩阵, 那么, 由矩阵理论可知, A (G) 的所有特征值均为实数, 设其特征值为λ1 (G) ≤λ2 (G) ≤…≤λn (G) .设D (G) 表示G的顶点的度对角矩阵, L (G) 表示G的拉普拉斯矩阵, 显然L (G) 也是一个实对称矩阵, L (G) 的最大特征值称为图的拉普拉斯谱半径, 记为μ (G) .根据Gerschgorin定理, 我们知道L (G) 的特征值是非负实数;又由它的行和为0, 可知0是L (G) 的最小特征值, 因此我们不妨设L (G) 的特征值0=μ1≤μ2≤…≤μn, 这里, 谱半径μ=μn.

二、循环图的背景

设S是集合{1, …, n-1}的任意一个子集, 并且满足S=-S mod n, 设图G具有顶点集{0, 1, …, n-1}, 若它的任意两个顶点i与j相邻当且仅当i-jmod n∈S, 则称图为循环图.显然循环图的邻接矩阵是一个循环矩阵, 即ai, j=ai-1, j-1.

一般地说, 循环矩阵具有下列形式:

$C\{\begin{array}{cccccc}c{}_0& c{}_1& c{}_2& \cdots & c{}_{n-2}& c{}_{n-1}\\ c{}_{n-1}& c{}_0& c{}_1& \cdots & c{}_{n-3}& c{}_{n-2}\\ & & \cdots \cdots & & & \\ c{}_2& c{}_3& c{}_4& \cdots & c{}_0& c{}_1\\ c{}_1& c{}_2& c{}_3& \cdots & c{}_{n-1}& c{}_0\end{array}$

由于循环图的邻接矩阵是一个实对称的循环矩阵, 而C的特征值为λj=${\displaystyle \sum _{k=0}^{n-1}}$$c{}_{k}e{}^{-\frac{2\text{π}jk}{n}i}‚j=1,2,\cdots ,n‚i{}^2=-1①.$

三、一些引理

引理1 若G是一个r正则图, 则μi-r=-λn-i+1, i=1, 2, …, n.

引理2 若G是一个具有n个顶点的连通图G (不是完全图) , 则$\frac{\mu {}_n}{\mu {}_2}\ge \frac{\text{Δ}(G)+1}{\delta (G)}$.由于循环图是一个正则图, 因而不等式$\mu {}_n\ge \frac{r+1}{r}\mu {}_2②$成立.

引理3 设G是一个含有n个顶点, m条边, 直径为D的简单图, 则$\mu {}_2(G)\ge \frac{2n}{2+n(n-1)D-2mD}③.$

等式成立当且仅当G=P3或G是一个完全图.

四、主要结论

定理1 若G是一个顶点为n, 度为r的循环图, 则

$\mu {}_n\ge \left(1-\underset{1\le k\le n-1}{{\displaystyle \max }}\left\{\cos \frac{2\text{π}k}{n}\right\}\right).$

证明 由于G是一个循环图, 那么根据1式, G的特征值为λj=${\displaystyle \sum _{k=0}^{n-1}}$$c{}_{k}e{}^{-\frac{2\text{π}jk}{n}i}.$

由引理1, L (G) 的特征值, 得μn=r-${\displaystyle \sum _{k=0}^{n-1}}$$c{}_k\cos \frac{2\text{π}k}{n}.$

令$s=c{}_0+c{}_1\cos \frac{2\text{π}}{n}+\cdots +c{}_{n-1}\cos \frac{2(n-1)\text{π}}{n}$, 由于c0, c1, c2, …, cn-1中ci (i=0, 1, …, n-1) 只能是0和1, 并且只有r个1, 其余都是0.

因此有$s=c{}_0+c{}_1\cos \frac{2\text{π}}{n}+\cdots +c{}_{n-1}\cos \frac{2(n-1)\text{π}}{n}\le$

r$\underset{1\le k\le n-1}{{\displaystyle \max }}$$\left\{\cos \frac{2\text{π}k}{n}\right\}$, 所以, $\mu {}_n\ge r-s=r\left(1-\underset{1\le k\le n-1}{{\displaystyle \max }}\left\{\cos \frac{2\text{π}k}{n}\right\}\right).$

定理2 若G是一个顶点为n, 度为r的循环图, 则

$\mu {}_n\ge \left(1+\frac{1}{r}\right)\frac{2n}{2+n(n-1)D-nrD}.$

证明 因为G是一个顶点为n, 度为r的正则循环图, 所以, $\mu {}_n\ge \left(1+\frac{1}{r}\right)\frac{2n}{2+n(n-1)D-nrD}④$.下面分三种情形予以讨论:

(1) 若D=1, 由于G是循环图, 所以G只能是完全图, 从而r=n-1, 因而由③式有, μn=n.

(2) 若D=2, 不可能是完全图, 因此, 由④式, 得

$\mu {}_n\ge \left(1+\frac{1}{r}\right)\frac{2n}{1+n(n-r-1)}.$

(3) 若$D\le \left[\frac{{\rm N}}{R}\right]\le \frac{{\rm N}}{R}+1$, 根据④式, 得

$\begin{array}{l}\mu {}_n\ge \left(1+\frac{1}{r}\right)\frac{2n}{2+n(n-r-1)D}\ge \\ \left(1+\frac{1}{r}\right)\frac{2nr}{2r+n(n+r)(n-r-1)}\ge \frac{2n(r+1)}{2r+\left(n-\frac{1}{3}\right){}^3}.\\ \end{array}$

参考文献

[1]B.Mohar, Some applications of Laplace eigenvalues of graphs[J].NATO ASI Series C, Vol.497, Kluwr, 1997.

[2]D.Cvetkovic, M.Doob, H.Sachs, Spectra of graphs:Theory and applications[M].3rd revised and enlargededition, J.A.Bart verglas, Heidelberg, Leipzig, 1995.

Laplace 篇4

Laplace算子是最典型的椭圆形偏微分算子,它在学习多元微积分部分就出现了,在直角坐标系它的形式非常简单,但是在极坐标下却出现了比较复杂的形式.如何将直角坐标的形式转化为极坐标形式是多元微积分的一道非常好的练习.我们将给出一个相对比较简单的一种几何的方法.这种方法是将Laplace算子看成是一种特殊Riemann度量下的Laplace算子.

1 Riemann度量下的Laplace算子

我们首先来引入Riemann度量下的Laplace算子的定义,这是一个一般性的定义,对任何坐标系都成立,因此十分重要.

定义1.1[1]设在任何坐标(x1,x2,…,xn)下,若有Riemann度量

其中G=(gij)构成一个矩阵,则相应的Laplace算子为

其中,g=det(gij),gij=(gij)-1.

2 直角坐标系下的Laplace算子

定义2.1[2]在直角坐标下(x,y,z)下,Laplace算子的形式为

如果用Riemann度量下的Laplace算子的定义,它实际上是如下度量的Laplace算子

相当于

我们将gij=δij代入方程(2)中,这时g=1,gij=δij,通过计算我们得到

3 极坐标下的Laplace算子

如何将直角坐标形式的Laplace算子(3)转化为极坐标形式是多元微积分的一道练习,在多元微积分的时候我们只能将极坐标用直角坐标表示,然后通过复合函数的求导得出极坐标的形式,虽然想法比较简单,但是具体的计算非常复杂.这里我们通过给出度量(4)的极坐标形式,然后就可以计算出相应的Laplace算子.

令

求微分得

代入到Riemann度量(4)中,得到

相应于该度量的Ｌａｐｌａｃｅ算子为

这样我们就得到了Laplace算子在极坐标下的表示.

摘要：Laplace算子在直角坐标系下有着非常简单的形式,但是在极坐标下却比较复杂,我们给出从直角坐标系下的Laplace算子到极坐标形式变换的一种几何方法.

关键词：Laplace算子,Riemann度量,极坐标

参考文献

[1]丘成桐,孙理察.微分几何[M].北京:科学出版社,1988.

Laplace 篇5

伴随计算机网络技术的迅速发展,数据传送通道的公共、信息发布平台的公开、数据高度的共享、以及数字化数据的便捷复制,极大提高了空间信息服务的质量,同时亦增加了数据的不安全隐患。航空航天技术的快速发展,促使遥感影像数据逐渐成为了民用空间服务的主要信息来源。数字化遥感影像的版权归属、数据的完整性已经成为可信赖空间信息服务的关键性技术问题之一。数字水印技术作为版权保护的有效技术手段[1],在图像、音频及视频等领域已经取得了较好的研究和应用,因此,遥感影像的数字水印技术也成为了实现遥感影像版权保护的可行的解决方案。

制约数字水印技术发展的最大技术难题是水印攻击[2],数字水印就是在不断的“攻”与“防”中发展起来的,目前为止,能够真正经受得住任何攻击的算法还不存在。在各类水印攻击中,几何攻击极为易行、难以抵抗和防御,且严重影响数字水印的有效性。因此,近年来针对数字水印的抗几何攻击鲁棒性,学者们陆续提出了大量的方法。文献[3]提出将小波应用于遥感影像特征的自适应二维盲水印算法,将水印灰度图像实施置乱加密和小波压缩之后,利用相邻特征的平均值和奇偶判决法嵌入遥感影像的二阶小波变换域内的所选子带上,水印检测无需原始遥感影像。文献[4]提出基于内容的离散余弦变换域自适应遥感图像数字水印算法,将二值水印图像安全嵌入遥感图像(纹理区)的离散余弦变换域内的高频系数中,对信号处理、几何剪切具有较好的抗差性。文献[5]提出了DFT域数字水印算法在遥感图像中的应用,该算法对水印灰度图像实施置乱后,采用一种基于内容自适应的离散傅立叶变换域的水印嵌入算法,将数字水印嵌入遥感图像的纹理区中以实现版权保护,并同时保证边缘信息不受较大损害。

上述文献分别提出了DWT、DCT和DFT域的遥感图像水印算法,算法各具优点,但均为基于图像全局的水印算法,抵抗行列剪切、任意角度的旋转、水平镜像等复杂几何攻击能力较差,仅要几种简单攻击组合就会彻底破坏水印的鲁棒性。

本文在上述研究的基础上,结合“第二代数字水印”思想,在遥感影像中选取特征区域冗余嵌入数字水印,提出了一种基于Harris-Laplace特征区域的遥感图像水印算法。

1 特征区域的确定

1.1 局部化水印基本思想

针对全局水印算法存在的缺欠,Kutter等提出了 “第二代数字水印”概念,基本思想如下:在图像中选择相对稳定的特征点区域嵌入水印,水印被独立嵌入每一区域内,水印定位与检测也利用特征点实现,实现了对几何攻击有效抵抗。第二代数字水印是基于图像内容的局部化水印技术,一般采用将数字水印冗余地嵌入图像中多个局部区域的方法,通过某一个局部区域就能够检测到水印的存在。

此类局部水印算法的稳定性依赖于:选取特征空间的方法(点特征、线特征、以及区域特征)、相似性度量的算法、搜索空间和搜索策略等[6,7,8]。选取到相对稳定的特征点是实现匹配校正、标识水印嵌入位置,以及实现水印定位与检测的重要基础。

1.2 Harris-Laplace特征点检测

Harris算子是一种有效的角点检测算子,在其基础上改进得到Harris-Laplace算子,Harris算子特点如下:Harris算子的计算仅用灰度一阶差分,易于实现;Harris算子通过计算图像每个点的兴趣值后,在邻域内选择最优点,因此提取的点特征均匀合理;Harris算子提取的点特征能够有效抵抗图像旋转、灰度值变化、噪声以及视点变换等。

这些良好特性促使Harris特征点在机器视觉、模式识别等方面得到广泛应用。然而,对于尺度变化较大的视觉系统,Harris保持特征不变性的能力较差,即对图像缩放变换的重复率低,使其对纵横比变换等攻击的鲁棒性较差。

特征点的尺度信息通常是未知的,一般而言,较大尺度能较可靠地消除误检,但定位不够准确;反之,较小尺度能准确定位真正的特征点,而误检率较高。因此,在Harris角点检测的基础改进得到Harris-Laplace算子,通过引入尺度空间理论,在多尺度空间上寻求最大响应点。具体实现过程如下:

(1) Harris角点的尺度空间表示如下:

${\rm M}(x,y,\delta {}_{{\rm I}},\delta {}_D)=\delta {}_D^2.G(\delta {}_{{\rm I}})*\left[\begin{array}{cc}L{}_x^2(x,y,\delta {}_D)& L{}_{x}L{}_y(x,y,\delta {}_D)\\ L{}_{x}L{}_y(x,y,\delta {}_D)& L{}_y^2(x,y,\delta {}_D)\end{array}\right](1)$

其中, (x,y)为图像中像素点的坐标,δI、δD分别为积分尺度和微分尺度,高斯尺度空间用La表示,计算函数在方向的偏导数,若δD已知,可定义如下:

La(x,y,δD)=Ga(x,y,δD)*I (2)

其中,数字图像记为,I均值为0、方差为δD时的高斯函数用G表示,*表示线性卷积。依据给定的δI与δD,来确定点(x,y)的梯度因子:

R(x,y,δI,δD)=Det(M(x,y,δI,δD))-η.Tr2(M((x,y,δI,δD))) (3)

其中,Det(·)为矩阵行列式,Tr(·)为矩阵的迹,η是常数(通常取值0. 04)。

条件1为:$R(x,y,\delta {}_{{\rm I}},\delta {}_D)>R(\widehat{x},\widehat{y},\delta {}_{{\rm I}},\delta {}_D)\forall (\widehat{x},\widehat{y})\in Q$

条件2为: R(x,y,δI,δD)≥tu

其中, Q为以点(x,y)为中心的邻域,tu是阈值。当点(x,y)的梯度因子R(x,y,δI,δD)同时符合上述条件时,表示该点在正交方向上梯度变化显著,将其作为候选特征点。

(2) 确定候选特征点集之后,计算每一特征点的特征尺度。特征尺度(characteristic scale)反映局部图像特征与操作算子间的最大相似度,即求特定尺度范围内,某函数极值点对应的尺度,本文选用LOG函数,计算特征尺度。LOG算子如下:

$L{\rm O}G(x,y,\delta {}_D)=\delta {}_D^2\left|\frac{\partial {}^{2}G(x,y,\delta {}_D)}{\partial x{}^2}+\frac{\partial {}^{2}G(x,y,\delta {}_D)}{\partial y{}^2}\right|*{\rm I}(4)$

(3) Harris-Laplace算子的实现过程:

对于给定的尺度空间δ${}_{{\rm I}}^{(n)}$=1.4nδ0,δ${}_D^{(n)}$=0.7δ${}_{{\rm I}}^{(n)}$(n=1,2,…,15)及阈值tu=1000,先通过尺度自适应Harris算子,求得预选特征点集{pk},然后,对于集合中的每个预选特征点,迭代确定最终的图像特征点与特征尺度。实现方法如下:

Step1 检验LOG算子在预选的图像特征点处,能否在尺度搜索范围内获得局部极值,如果得不到极值,该点舍弃。尺度搜索范围限定为:δ${}_{{\rm I}}^{(k+1)}$=tδ${}_{{\rm I}}^{(k)}$,其中t=0.7,…,1.4。

Step2 将LOG算子取极值的特征点pk加入候选特征点集合{p′k}。

Step3 重复以上两步,至所有特征点处理完毕为止。集合{p′k}就是所求Harris-Laplace特征点。

基于Harris-Laplace特征点抗压缩、噪声、滤波、旋转、缩放等的良好特性,根据所选的特征点及其特征尺度,确定水印的嵌入区域,提高算法的抗几何攻击能力。

1.3 几何不变特征区域的确定

依据得到的特征点集建立几何不变区域作为水印嵌入位置,几何不变区域能够有效抵抗绝大多数几何攻击,特别是局部几何攻击。几何变换不变区域的选择遵循如下原则:几何不变区域可逆,含水印图像遭受攻击时,几何不变区域依然能被提取出,确保水印检测正确;几何不变区域的自适应,水印图像遭遇几何攻击,不变区域应能随之自适应地变换,使其中内容保持不变,正确检测水印[9]。

依据上述原则,本文确定几何不变区域方法如下:将水印嵌入的几何不变区域确定为以特征点为圆心,kσ为半径圆片,k取常数(6-10)。特征点的特征尺度用σ表示,σ过大、过小的特征点舍弃,若嵌入范围重合,则舍弃σ较小的特征点,将多尺度上的稳定性较好的特征点保留,进而确定几何不变区域。

2 基于Harris-Laplace特征区域的遥感图像水印算法

2.1 数字水印的嵌入

设原始载体遥感图像1大小为M×M,水印嵌入具体步骤如下:

(1) 水印生成

通过种子密钥K,生成一个随即序列,对其进行二值化,选取长度为N的一段,作为数字水印W。记为W={wi,1≤i≤N}。

依据1.2节方法,提取Harris_Laplace特征点,依据1.3节方法,并构建几何不变区域。

(2) 嵌入水印信息

依据待嵌水印位数N,将每个圆形区域划分成若干同心圆环,如图1所示。在每个所选尺度不变区域冗余嵌入水印。

首先,把像素坐标(x,y)变换为极坐标:极坐标半径ρ,角度θ。由圆心起的第i个圆环区域CRi记为:

$CR{}_i=\left\{(x,y)|(i-1)\cdot \frac{R}{{\rm N}}\le \rho {}_{x,y}<i\cdot \frac{R}{{\rm N}}\right\}$ (5)

其中 i=1,2,…,N。

然后,依据待嵌入的水印信息,采用奇偶量化方法,由圆心向外修改每个圆环的像素,将圆环修改为“奇圆环”或“偶圆环”。利用公式(6),对每一个像素I(x,y)赋值“0”或“1”,Δ为量化步长。

$Q(x,y)=\{\begin{array}{l}0k\text{Δ}\le {\rm I}(x,y)<(k+1)\text{Δ}k=0,\pm 2,\pm 4,\cdots \\ 1k\text{Δ}\le {\rm I}(x,y)<(k+1)\text{Δ}k=\pm 1,\pm 3,\pm 5\cdots \end{array}(6)$

为增强水印鲁棒性,将修改后的像素值限定于相应量化区间中间值。为此,本文依如下规则对像素进行修改。

${\rm I}{}^{*}(x,y)=\{\begin{array}{l}r\cdot \text{Δ}+0.5\text{Δ}Q(x,y)=w{}_i\\ r+1.5\text{Δ}Q(x,y)\ne w{}_{i}r\cdot \text{Δ}>0.5\text{Δ}\\ r-1.5\text{Δ}Q(x,y)\ne w{}_{i}r\cdot \text{Δ}\le 0.5\text{Δ}\end{array}(7)$

其中$r=\lfloor \frac{{\rm I}(x,y)}{\text{Δ}}\rfloor$,在每一所选圆形区域冗余嵌入水印,得到含水印遥感影像,水印嵌入核心流程如图2所示。

2.2 数字水印的检测

水印提取是水印嵌入的逆过程,本文的水印提取无需原始遥感影像参与,实现步骤如下:

(1) 同步水印信息

提取Harris-Laplace特征点,同步含水印的区域,得到一系列圆形区域。

(2) 水印信息提取

在每个同步的圆形区域内进行圆环划分,其方式与嵌入过程同步,依据公式(6),计算每一圆环区域内像素的量化函数值, NUMi,0表示函数值为0的像素数目,函数值为1的像素数目记NUMi,1,依据式(8)提取水印信息:

然后,计算原始水印w与提取水印w′的归一化相关函数(NC):

${\rm N}C=\frac{{\displaystyle {\sum }_{i=1}^{{\rm N}}w}(i)w^{\prime }(i)}{\sqrt{{\displaystyle {\sum }_{i=1}^{{\rm N}}w}{}^2(i)}\sqrt{{\displaystyle {\sum }_{i=1}^{{\rm N}}{w}^{\prime }}{}^2(i)}}$ (9)

每一圆形区域内均可提取到水印信息,选取其中NC值最大的作为最终的水印检测结果,水印提取核心流程如图3所示。

3 实验与结果分析

为了验证所提算法有效性、测试分析水印透明性、相关性、及抗攻击的鲁棒性,分别进行了特征点提取、水印的嵌入与提取实验、水印的常规信号处理和几何攻击等攻击实验。实验中,采用的原始宿主图像为截选256×256大小的遥感影像图,水印为通过种子密钥K,生成一个随即序列,对其进行二值化,选取长度为N的一段,作为数字水印W,水印长度64bit。按本文提出的算法进行了水印的嵌入与提取、不同参数的水印攻击及性能检测仿真实验。

3.1 Harris Laplace特征点提取

按照本文1.2和1.3节方法进行Harris Laplace特征点提取和特征区域确定。图4是宿主遥感影像提取的Harris-Laplace特征点以及构建的尺度不变区域:(a) 提取的Harris-Laplace特征点;(b) 所有的尺度不变区域;(c) 选定的不变区域。

Harris Laplace特征点提取核心编程实现:

3.2 水印性能检测

3.2.1 透明性检测

图5为嵌入水印前后的遥感影像,可见原始遥感影像与含水印影像视觉差别不明显,含水印影像与原始遥感影像峰值信噪比(PSNR)为41.24dB。

3.2.2 相关性检测

通过100个不同水印进行实验检测,以验证水印相关性,其中被嵌入的水印是第50个。仿真实验结果如图6。响应曲线仅在第50的位置有明显峰值,可见检测器对真实水印的响应远高于其它错误水印。

3.2.3 鲁棒性检测

对含水印遥感图像进行常见的几何攻击包括:旋转、缩放、平移、剪切,纵横比改变等。对于平移、剪切攻击,只要保留一个几何不变区域,水印即可正确提取。平移、剪切攻击结果示例如图7,提取水印NC值均等于1;尺度缩放、旋转攻击如图8、图9,图为将水印分别嵌入两幅遥感影像中的水印检测NC值与缩放因子、旋转角度的关系。

3.3 实验分析

通过上述实验可知,本文算法对RST、非等比例缩放、局部剪切等常见几何攻击具有较好鲁棒性。

(1) 基于在图像每一局部区域嵌入水印,对平移、剪切攻击,只要保留任一局部区域即可成功检测水印;

(2) 鉴于水印嵌入的区域具有缩放不变性,进行缩放攻击时,局部区域内的内容是相同的,采用统计提取水印位的方法,即使像素有所变动,只要该区域内一半像素改变不大,亦可正确提取水印位;

(3) 由于Harris-Laplace特征点对旋转有良好稳健性,且本文选取圆环区域嵌入水印,因此具有较好的旋转不变性。

4 结 语

现有DWT、DCT和DFT域的遥感图像水印算法,虽各具优点,但均为基于图像全局的水印算法,抵抗行列剪切、任意角度的旋转、水平镜像等复杂几何攻击能力较差,仅要几种简单攻击组合就会彻底破坏水印的鲁棒性[10,11]。本文针对上述问题,结合“第二代数字水印”思想,在遥感影像中选取特征区域冗余嵌入水印信息。算法通过种子密钥K,生成一个随即序列,对其进行二值化,选取长度为N的一段,作为数字水印信息;然后,在宿主遥感影像中提取Harris-Laplace特征点,并依据其确定水印嵌入特征区域;最后根据待嵌入的水印信息,将所选特征区域划分为若干圆环,采用奇偶量化方法嵌入数字水印。实验表明本文算法对JEPG压缩、噪声与滤波等常规信号处理具有较好的鲁棒性,同时抵抗旋转、缩放、平移及其联合等几何攻击。

参考文献

[1]陈丽萍,姚志强,孔祥增.基于证据理论数据融合的双效水印[J].计算机应用与软件,2011,28(8):50-53,161.

[2]Bami M,Cox I J,Kalker T.Digital watermarking[C]//Lecture Notes in Computer Science(LNCS),USA(New York):Springer-Verlag,2005:15-17.

[3]王贤敏,关泽群,吴沉寒.小波用于基于遥感影像特征的自适应二维盲水印算法[J].计算机工程与应用,2004,(20):37-41.

[4]王向阳,杨红颖,邬俊.基于内容的离散余弦变换域自适应遥感图像数字水印算法[J].测绘学报,2005,34(04):324-330.

[5]陈晨,何建农.DFT域数字水印算法在遥感图像中的应用[J].计算机与现代化,2008,159(11):77-79.

[6]Dong P,Galatsanos N P.AffineTransformation Resistant Watermarking Based on Image Normalization[J].IEEE Trans.on image processing,2002,8:865-882.

[7]Li Leida,Guo Baolong,Shao Kai.Geometrically robust image water-marking using scale-invariant feature transform and Zernike moments[J].Chinese Optics Letters,2007,5(6):332-335.

[8]孟岚,杨红颖,王向阳.基于图像归一化的DFT域数字水印嵌入算法[J].小型微型计算机系统,2008,29(11):2153-2156.

[9]肖迎杰,孙劲光,张云.一种模板和图像特征相结合的数字水印算法[J].计算机应用,2010,30(1):50-53,57.

[10]Lee S,Choi H,Choi K,et al.Fingerprint-quality index using gradient components[J].IEEE Transactions on Information Forensics and Secu-rity,2008,3(4):791-799.

Laplace 篇6

方程组问题及变分问题的研究, 一直是一个广泛的研究课题, 针对不同的情况, 很多数学工作者利用不同的方法研究了该类问题, 做出了许多有意义的工作 (见文[1-3]及其参考文献) , 它的应用背景见[4,5]。

本文的主要工作是采用锥不动点定理研究了在单位球 上的 方程组:

的正径向解的存在性, 其中, ai (x) 和p (x) 是可微的径向对称函数

本文的模型 (1) 推广了文[1]中的模型 (1) , 并减弱了文[1]中的条件λ→0。

2 准备工作

方程组 (1) 的径向解满足

其中, 且满足方程组 (2) 。

定义1若向量值函数中, 则称解是正的, 本文将证明方程组 (2) 至少有一正的非平凡解。

引理1[6]E是赋予范数||·||的Banach空间, K是E中的锥, 对r>0, 定义为Kr的对于K的边界, 假设是完全连续的。

(I) 若存在使得则有

(II) 若, 则有。

3 主要定理及引理

定理1假设是连续的, 且满足, 则方程组 (1) 有一个正的径向解。

设为Banach空间, , 定义范数

定义为X中的锥, 设是一组的映射组, 定义:

引理2是紧的和连续的。

证明 (1) 由于是连续的且, 再由T的定义可得。

(2) 对K中有界集D, 存在常数r, 使得对于任意, 有.由于fi (u) 是连续的, 则存在d, 使得, 记设是D中的一有界数列, 则有

由此得T (D) 是一致有界的。

对使得对每一有

所以, T (D) 是等度连续的。由Arzela-Ascoli定理得T (M) 是列紧的, 因此是紧的。

(3) 设{um}是K中任意收敛于的序列, 且有使得对所有, 有, 假设Ti不连续, 则和{um}的一子列{umj}使得

在[0, 1]上连续, 在闭集上连续, 由控制收敛定理有

由Ti的紧性得一{umj}子列 (仍记为{umj}) , 使得Tiumj一致收敛于, 因此由 (3) 式得

由逐点收敛性 (4) 式可得, 这与上式 (5) 矛盾, 因此Ti连续, 证毕。

4 定理1的证明

由引理2可得:T是完全连续的.

若, 则有

故由上式 (6) 选取适当的r2>0, 使得

由已知条件, 则和常数η使得

其中且, 选择η使得, 时满足和时满足

若, 则得方程组 (1) 的解, 因此本文考虑对所有, 取于于于, 因此。本文证明, 对所有的都成立, 假设及使得, 下面找出矛盾:

因为, 由假设得t0>0, 又因为, 由v的定义, 得到对所有, 有, 设

从而得到且

现在对有

由范数的定义, 由式 (7) 对有

因此, 这矛盾于t*的定义, 从而由引理1得

由Leray-Schauder度的可加性得。因此T在中有一个不动点, 即为方程组 (1) 的正解, 定理1证毕。

参考文献

[1]任淑青, 刘辉昭.一类方程组的径向解的存在性[J], 河北工业大学学报, 2008, 37 (2) :56-59.

[2]Wang H.An existence theorem for quasilin-ear systems[J].Proceedings of the Edinburgh Mathematical Society, 2006, 49:505-511.

[3]X.L.Fan, X.Y.Han.Existence and multi-plicity of solutions for equations in[J].Nonlinear Analysis59 (2004) :173-188.

[4]M.ruzicka, Electrorheologi-cal Fluids:Model ing and Mathematical Theory[M], Spring, Berlin, 2000.

[5]V.V.Zhikov, Averaging of functionals of the calculus of vatiations and elasticity theory[J], Math.USSR.Izv.9 (1987) 33-66.

Laplace 篇7

针对以上问题, 本文提出了一种基于Harris-Laplace特征点检测的强韧性数字水印算法。

1 图像的Harris-Laplace特征点检测

Harris检测法是在图像处理中应用最广的一种方法[4,6]。它对于抵抗图像的旋转和平移都有高度的抵抗力, 但不具有尺度不变性, 所以, Mikolajczyk[6]等学者则针对这个问题提出改善方法, 借由建立尺度空间 ( scale-space) 来侦测角点更能准确定位出真正的角点, 此方法是利用M矩阵估计特征的强度, 通过阈值则视为角点。M矩阵的定义如式 ( 1) 。

式 ( 1) 中, G ( σI) 表示高斯函数, σI为积分尺度, σD微分尺度, x, y为图像坐标, L ( x, y) 则表示计算高斯平滑影像x与y方向导数。

由上述公式可知, 对于给定的 σI和 σD, 可以确定出像素点x, y的特征强度, 如公式 ( 2) :

式 ( 2) 中R代表角点强度、det表示矩阵的行列式、trace代表矩阵的迹, k为常数, 一般设为0. 04。当像素点的特征强度大于阈值, 且该像素点为尺度空间的极值时, 就可以认为该像素点为图像的角点。

对于各尺度的空间极值点, 该极值点的尺度空间则代表图像的局部特征, 可与其他尺度空间极值进行运算, 进而萃取出具有尺度不变特征的局部特征点。如何萃取出具有尺度不变的特征点则使用LOG ( Laplaceian-of-Gaussian ) 来做运算, 而LOG运算是将高斯滤波和拉普拉斯边缘检测做结合, 定义如公式 ( 3) 所示。

只要给定图像的像素点与尺度范围, LOG运算局部极值所对应的尺度即为特征尺度, 就得到具有尺度不变的特征点。具体步骤如下:

( 1) 设定好尺度空间范围。

σ (n) I=1.2n×1.5, σ (n) D=0.7σ (n) I (n=1, 2, …, 13) 与阈值 (T) , 利用公式 (1) 和 (2) 计算出每一尺度空间上的候选特征点Pk。

( 2) 使用迭代法验证每个尺度空间候选特征点Pk的LOG运算是否在整个尺度空间搜索范围内获得区域极值, 如不能获得极值, 则舍弃该点, 搜寻范围限定在 σI (k+1) = tσI (k) , 其中t=0. 7, 0. 8, …, 1. 4。

( 3) 对于LOG运算能获得极值图像特征点pk, 在该点区域内计算特征强度R最大的特征点pk+1, 若pk+1存在则舍弃pk。重复b, c过程直到 σI (k+1=σI (k) 或者pk+1不再变化。

2 数字水印的嵌入

选用的嵌入法人水印图像是有意义水印, 首先要对水印图像进行预处理, 将有意义水印图像经过Arnold置乱转变为水印序列后利用重复嵌入的方式将水印序列嵌入到原始图像中。水印图像经过处理后, 提取原始图像的特征点并计算质心与最显著特征点, 通过特征点的重定位来确定嵌入位置, 修改特征点与质心和最显著特征点的角度以嵌入水印序列信息, 之后通过德劳内三角技术对修改后的特征点进行调整, 利用仿射变换将改变后的图像矩阵转换为原先设定的矩阵, 最后完成水印信息的嵌入过程, 得到含水印的图像信息。

该嵌入过程的示意图如图1 所示。

具体嵌入步骤如下:

步骤1 用Harris-Laplace角点检测方法提取特征点, 先建立一个尺度空间, 利用Harris角点检测出图像中区域的极值点作为候选点, 最后利用迭代法计算所有尺度范围内候选点的LOG ( Laplacian-ofGaussian) 得到特征点的位置和尺度空间。

Harris-Laplace的角点检测提取的特征点信息包括x和y坐标位置、尺度空间 ( σ) 与特征强度 ( R) 。因图像的特征点是在不同的尺度范围上, 一般认为在较大的尺度下能找出真正的特征点, 但特征点的位置较不易确定。相反来说, 较小的尺度对特征点的位置较准确, 但容易找出错误的特征点, 所以, 本文将各特征点的尺度 ( σ) 、特征强度 ( r) 和质心到特征点的距离d相除并做降幂排序, 得到既稳固又强韧的特征点。如图2 所示。

但提取的特征点可能相邻太近或是靠近边界造成之后取出的特征点错误, 如此, 可借由x和y坐标与尺度 ( σ) 就能把太相邻或靠近图像边界的特征点删除。首先, 将尺度空间做降幂排序, 以固定倍数的尺度 ( σ) 半径之点坐标为中心点画图, 计算出外接正方形四个顶点。利用四个顶点x和y坐标, 来判断与其他特征点坐标是否有重叠。如图3 所示。

步骤2 计算质心和最显著特征点, 修改这些特征点和质心、最显著特征点与特征区域平均值, 来嵌入水印信息。

虽然取得一定数量的特征点, 但图像遭到攻击后特征点会变动。因此本方法找出定位整张图像特征点质心o与最显著特征点p, 质心o则是利用图像的质量重心 ( center of mass) 。其计算方式如公式 ( 4) 。令图像f ( x, y) 的大小为Mx N, 其 ( p+q) 阶的二阶矩阵 ( moment) :

式 ( 4) 中心矩 ( central moment) 定义为:

式 ( 5) , ( x0, y0) 即是图像f ( x, y) 的质心为o, 此点也具有几何不变的特性。

另一个则是最显著特征点p, 方法是将各特征点的尺度与其特征强度相乘再除以质心至特征点距离d, 得到的最大值的特征点来定位整张图像的位置。

式中 σi、Ri为i个特征点的尺度与特征强度索引, di为质心到特征点第i个距离。

质心o与最显著特征点p除了定位整张图像外, 这两点也是划分数字水印嵌入的区域, 连接质心o与最显著特征点p两点, 每360° / m顺时针划分成m个区域, 每个区域嵌入一个单位, 数字水印W是经过转换后的0, 1 序列, 数字水印长度为划分的m个区域, 位于每个区域的特征点均嵌入相同的单位。

步骤3 重新定位特征点, 修改特征点与质心、最显著特征点形成的角度及特征点尺度 ( σ) 大小的区域平均值来进行嵌入, 角度的嵌入过程如下:

设质心o、最显著特征点p与特征点q、计算出质心o与最显著特征点p向量记为; 将质心o分别和特征点q连接其向量记为, 使用cos-1函数可得其夹角 ( θ) , 如公式 ( 6) 。

特征点所求得的角度则使用一个嵌入函数进行水印的嵌入, 如公式 ( 7) 所示。首先将要嵌入的水印图像利用转换函数转换为0, 1 的水印序列来进行嵌入。假设我们所要嵌入的水印为0 的话, 则将特征点构成的夹角变为4 的倍数; 反之嵌入的水印为1, 则将角度改为4 的倍数加2. 依序将特征点进行计算可得到嵌入后的特征点的新角度

式 ( 7) 中 θi为特征点与质心、最显著特征点的夹角, i为第i个特征点, α 为水印的权重。α 值越大嵌入的水印越强韧, 但对图像的影响也越大。

由于是通过改变角度进行的嵌入, 其坐标也会跟着移动, 为避免特征点坐标改变太大, 造成特征点位置错误和图像失真严重, 要在一个n×n的范围内寻找改变角度后的特征点的新坐标来避免特征点严重移位。举例来说, 假设要嵌入的特征点为q, 其坐标为 ( 273, 115) 与质心o、最显著特征点p形成的夹角 θ 为53 度, 要嵌入的水印为1, 即嵌入水印后其夹角应为54°, 为避免特征点移位太多, 本文将特征点移动范围控制在坐标周围的2×2 区域内, 找出最接近角度的特征点坐标, 得到嵌入水印后的新特征点q。如图5 所示。

特征区域平均值嵌入方法如下:

根据Harris-Laplace detector获得尺度 ( σ) , 特征区域半径r的大小由公式 ( 8) 定义, round表示四舍五入, τ 为正整数的倍数, 这里定为3 的缩放倍数。

在求得特征区域后, 嵌入方法定义如公式 ( 9) 所示。假设特征点角度嵌入的水印信息为1, 其特征区域平均值也变为4 的倍数加2。通过角度与特征区域平均值的嵌入可以增加水印提取的正确性。

式 ( 9) 中为特征点与质心、最显著特征点的特征区域平均值, i为第i个特征点, α 则为水印的权重。

步骤4使用Delaunay三角网络技术, 将相邻的特征点每三个点形成一个三角形区域, 而全部特征点所构成的三角形样板区域T则包含各三角形特征点的顶点和顺序, 因微调后的特征点连接三角形顺序可能会与原始特征点所形成的三角样板不同, 只要利用原始特征点产生的三角形样板应用在位移特征点上, 就可保证原始与新特征点所形成的每个三角形网格区域是相同的, 如图6 所示。

步骤5 进行仿射转换。将更改坐标位置的特征点利用仿射转换将各特征点所形成的三角形区域转换至预先设定好的与原图大小相同的矩阵, 其余各点依序对应进行计算, 运算完后即完成了整个数字水印的嵌入过程。

3 数字水印的提取

数字水印的提取过程与嵌入过程相似, 其提取过程示意图如图7 所示。

首先, 利用先前的Harris-Laplace提取特征点, 再分别找出其质心与最显著特征点, 求其余特征点与质心和最显著特征点的角度与特征点区域平均值, 并利用提取函数进行水印提取, 其公式如下:

式 ( 11) 中wai为取出的角度水印序列, wvi为取出的区域平均值水印序列。

然后, 根据公式 ( 12) 计算各区域特征点角度与特征区域平均值为1 或0 的个数, 即可取出嵌入的水印信息Wi。

式 ( 12) 中wi为提取的水印特征点, 其中i表示为第i个特征点。

4 实验结果分析

本实验使用512×512 的Lena图像为原始图像, 如图8 ( a) 所示, 水印图像为64×64 的标有天津财经字样的二值图像, 如图8 ( b) 所示, 并且利用Matlab[6]进行仿真实验。以信噪比PSNR ( peak signal to noise ration) 来作为原始图像与嵌入水印图像质量的评估标准, PNSR值越大, 则图像差异度越小。嵌入后的图像如图8 ( c) , 显然水印的不可见性较好, 其PNSR值为31. 26d B, 图8 ( d) 是提取的水印, 它与原始水印的归一化相关系数NC值为1。表1 是本算法对常见图像处理和几何攻击的抵抗能力。并与Lee等学者[4]的方法作比较。

由实验结果可看出, 不管遭受何种攻击, 只要能正确取出质心与缩放比例乘以特征强度最大的点, 数字水印就很容易嵌入, 并且通过比较可知, 本文算法能很好地抵抗压缩、旋转、缩放等攻击, 具有较好的不可见性和鲁棒性。

5 结束语

本文利用Harris-Laplace提取的特征点与质心和Delaunay三角网络技术方法将特征点构成三角形, 将数字水印嵌入其中, 再利用正规化的概念经由仿射转换后得到嵌入数字水印的图像。由实验结果可知, 本研究是一个不需要原图就可将水印取出的方法, 且图像遭到几何攻击后也有高度的抵抗能力, 对于一般图像处理也有显著的效果。

参考文献

[1] Lauzef L M.The improbability of Harris interest points.IEEE Transactions on Pattern Analysis and Machine Intelligence, 2010;32 (6) :1141—1147

[2] Wang X Y, Yang Y P, Yang H Y.Invariant image watermarking using multi-scale Harris detector and wavelet moments.Computers and Electrical Engineering, 2010;36 (1) :31—44

[3] 袁大洋, 肖俊, 王颖.数字图像水印算法抗击和攻击鲁棒性研究.电子与信息学报, 2008;30 (5) :1251—1256

[4] 柏均, 张敏瑞, 梁文莉.基于特征点的图像水印嵌入方案.计算机工程, 2011;37 (13) :139—140

[5] 赵洁, 武斌, 张艳.基于Harris兴趣点和空域均值信息的图像复制粘贴窜改检测算法.计算机应用研究, 2013;30 (9) :2971—2973

[6] Mikolajczyk K, Schmid C.Scale and affine invariant interest point detectors.International Journal of Computer, 2001;60 (1) :147—151

[7] 张强, 王正林.精通Matlab图像处理.北京:电子工业出版社, 2009