有点编程基础的童鞋看到这个标题可能会有点懵逼,这还是个问题吗?不就是个等号(=)解决问题嘛!我也希望是如此简单,因为上个星期被这个问题折磨到崩溃!
一般的python程序需要赋值时的确是通过等号(=)实现的,不管是变量还是数组,例如:
i=1
pi=3.1415926
x=numpy.arange(1,10)
也可以实现一些稍微复杂的操作:
ilon=90
jlat=40
corr = Corr_piont_and_plane( djf[:,jlat,ilon] , djf )
其中Corr_piont_and_plane是一个函数,有两个参数,djf是冬季距平(500hPa位势高度场),三维数组,函数实现的功能是选择其中一个点(参数1,一维数组)和这个冬季距平(参数2)求时间相关,返回一个二维的相关系数矩阵,表征相关系数的空间分布。研究气候的童鞋会经常用到这个函数,只不过参数1大概率是某个指数序列,参数2大概率是前冬海温距平。
如果这个Corr_piont_and_plane函数是python写的,以上的代码是没有问题的,返回的结果也是正确的。但如果玩了一点骚操作,Corr_piont_and_plane是由Fortran或者C写的,编译成python可调用的动态库(如果对这个骚操作感兴趣,可以点这里),那前面的代码就是个大坑,我在这个大坑里转了一个星期。我先直接给明确结论,这个代码应该怎么写,然后再解释为什么。正确的代码应该这么写:
ilon=90
jlat=40
point=djf[:,jlat,ilon].copy()
corr = Corr_piont_and_plane( point, djf )
这个就涉及到今天这个技术贴的主题,numpy的赋值,更加python或者numpy的说法是拷贝。等号(=)实现的是浅拷贝,而numpy.copy()实现的是深度拷贝。尽量用简单的语言进行解释:在一个python程序,新建一个numpy的array之后,内核会分配一组内存用来保存数据,但是怎么进行管理呢?内核认识这片内存的地址,但是这个地址保存着数据,数组的名字,比如djf,指向这个数组的首地址。打个比方,一个小区,有好几栋楼,每栋楼好几层,然后一梯两户或者一梯三户,每个房子对着一户人家,小区物业怎么管理呢,正常情况下是几栋几零几(地址)住着谁家(数据),然后邻居之间关系好,会拥有下一户业主家的秘钥(下一个数据的地址,至于下一户是楼上楼下,还是同层别的住户,甚至于别的楼同门牌号的业主,和内存优先机制有关,列优先还是行优先,不同语言不一样,就是看小区物业自己咋规定)。回归主题,当进行索引或者切片的时候,返回的其实是地址,不是数据,内核会访问内存地址获取数据。但是当切片出来一个数组,比如这个point,等于(=)相当于给这些数据打上标记,叫做引用,数据的内存地址没有变化,但是用深度拷贝,相对于创建一个新数组保存新的数据,内存地址发生变化了。下图是一个简单示例,b是a的浅拷贝,所以b就是a,但是c是a的深度拷贝,如果c不是a。还用小区做例子进行解释,就是有一栋楼的某一个扇面的下水道坏了,要进行处理。小区物业记下了这一个扇面属于哪一栋楼的哪个哪个扇面,安排维修师傅更换下水道,这属于浅拷贝,因为住户还都住在这里,当维修人员发现楼房年久失修,更换下水道也不管事,只能异地拆迁的时候,属于深度拷贝,因为住户换地方住了。
再次回归主题,当所有的代码都是基于python编写,直接把切片后的数组当做参数进行传递没有任何问题,但是C-API接口进行传递时,就会存在问题,问题在于参数传递的是地址,到了Fortran编写的subroutine里,访问数组的首地址能够正确识别,但是后面的数据就错了,因为传递过来的地址Fortran不认识,可能原因这么几点,还需要进一步确认:1、通过接口传递切片数组时,抹掉了切片信息,子程序默认还是内存的下一个地址去读数据,更深层次的原因可能是python对变量进行了对象封装,叫做PyObject,2、python和C默认行优先,Fortran默认列优先,数据寻址找错了方向。基于这些原因,需要进行深度拷贝,重建参数数组的内存地址逻辑,再进行参数传递,才是正确的操作。
最后,对结论进行总结:python&numpy的赋值一般通过浅拷贝实现,即等于(=),也可以把切片数组作为参数进行传递,代码简洁,节省内存,当需要使用C-API接口传递参数,建议对参数使用深度拷贝赋值之后再进行参数传递。