PyTables

chimez

2021-03-01 11:23

PyTables: 分层数据集

PyTables 是基于 HDF5 格式构建的数据集，比 pandas 更底层，比 h5py 更高层。适合用在存储比 csv 表格更复杂更大的数据，同时又不想自己实现一些基本的数据操作的情况中，并且经过简单的处理就能用 pandas 进行进一步数据分析。

安装

   # 使用 conda
   conda install pytables
   # 或者 pip
   python3 -m pip install tables

测试

   import tables as tb

   tb.test()

基本使用

pytables 的基本逻辑是

创建一个 HDF5 文件

    f = tb.open_file("filename.h5", mode="w", title="file title", filters=tb.Filters(complevel=9))

    # 别忘了还要关闭它
    f.close()

在文件中创建任意层的数据集

    g = f.create_group("/parent/group", "group_name", title="group title")

在数据集中创建若干个数据表

    # 需要创建对数据的描述
    class data_description(tb.IsDescription):
        name = tb.StringCol(20)
        idn = tb.Int64Col()
        class sub_data(tb.IsDescription):
            name = tb.StringCol(20)
            id2 = tb.Float64Col(shape=(3,5))

    # 然后创建空的数据表       
    t = f.create_table(g, "table_name", data_description)

    # 在里面循环添加记录
    for i in range(10):
        t.row["name"] = "a"
        t.row["idn"] = i
        t.row["sub_data"] = (f'asd_{i}', np.random.rand(3, 5))
        t.row.append()

    # 最后刷新一下缓存，确保数据都写入到磁盘
    t.flush()

遇到的问题

在不支持 `flock` 的集群上运行

需要设置环境变量

   export HDF5_USE_FILE_LOCKING=FALSE

h5py

chimez

2021-02-28 22:11

h5py: 储存数据

h5py 是 HDF5 库的 python 封装，基本提供了所有对应的 C API，因此适合用在简单并底层的应用中，特别适合对 HDF5 文件有定制的应用中。

存

   with h5py.File('data.h5', 'w') as hf:
       hf.create_dataset("mat1", data=mat1, compression="gzip")
       hf.create_dataset("mat2", data=mat2, compression="gzip")

       g = hf.create_group("index")
       g.create_dataset("mat1", data=mat1, compression="gzip")

取

   with h5py.File('data.h5', 'r') as hf:
       mat1 = np.array(hf['mat1'])
       mat2 = np.array(hf['mat2'])

       g = hf['index']
       g_mat1 = np.array(g['mat1'])

OpenMP

chimez

2021-02-20 15:37

OpenMP 共享内存并行

OpenMP 是单机多线程共享内存并行工具。

OpenMP 使用 fork-join 模型进行并行，即主线程通过 fork 创建多个并行线程，之后再将子线程 join 回主线程。

OpenMP API 由三部分组成

编译器指令 #pragma omp ...
运行时库函数 #include <omp.h>
环境变量 OMP_NUM_THREADS 等

参考资料

llnl/openmp 的教程

编译器指令

编译器指令的格式为

`#pragma omp`	指令名	`[clause, ...]`	换行
所有OpenMP指令的开头	一系列可用的指令	可选的，从句顺序任意	换行不需要分号,可用用反斜线折行

   #pragma omp parallel default(shared) private(beta, pi)

parallel 指令

创建一系列线程运行接下来的结构块，并在结构块结束时合并所有线程

可用的从句

if (scalar_expression)
private (list) 指定变量是私有的
shared (list)
default (share | none)
firstprivate (list)
copyin (list)
num_threads (integer-expression)

for 指令

指定接下来的循环迭代并行运行。

schedule (type [,chunk]) 描述循环如何分配给线程
ordered 指定线程必须按照循环的顺序合并
private (list)
firstprivate (list)
lastprivate (list)
shared (list)
reduction (operator: list)
collapse (n)
nowait 循环结束后不等待未结束的线程

    #pragma omp parallel shared(a,b,c,chunk) private(i)
      {

      #pragma omp for schedule(dynamic,chunk) nowait
      for (i=0; i < N; i++)
        c[i] = a[i] + b[i];

      }  /* end of parallel section */

sections, section 指令

外层用 sections 指令，内部用 section 指令

    #pragma omp parallel shared(a,b,c,d) private(i)
      {

      #pragma omp sections nowait
        {

        #pragma omp section
        for (i=0; i < N; i++)
          c[i] = a[i] + b[i];

        #pragma omp section
        for (i=0; i < N; i++)
          d[i] = a[i] * b[i];

        }  /* end of sections */

      }  /* end of parallel section */

single 指令

指定块内部同时只有一个线程运行，用来处理以下线程不安全的操作

task 指令

master 指令

指定只有主线程执行块内的程序。

critical 指令

指定同时只有一个线程能运行内部的代码。

    #pragma omp parallel shared(x) 
      {

      #pragma omp critical 
      x = x + 1;

      }  /* end of parallel section */

barrier 指令

要求所有线程立即同步

taskwait 指令

atomic 指令

执行原子操作，相当于小的 critical 指令

flush 指令

要求所有缓冲的修改写回内存中

下面一些指令隐含了 flush 指令

barrier
parallel 的进和出
critical 的进和出
ordered 的进和出
for 的出
sections 的出
single 的出

ordered 指令

要求代码块内部的指令顺序执行，只能用在 for 指令内部

threadprivate 指令

指定一些变量是线程私有的，每个线程初次使用这些变量时都是未初始化的

reduction 从句

要求循环分成相等的部分，并在全部运行最后再执行指定的指令

支持的格式

    x = x op expr 
    x = expr op x 
    x binop = expr
    x++ 
    ++x
    x--
    --x

    op: +, *, -, /, &, ^, &&, ||
    binop: +, *, -, /, &, ^, |

例子

    #pragma omp parallel for      \  
      default(shared) private(i)  \  
      schedule(static,chunk)      \  
      reduction(+:result)  

      for (i=0; i < n; i++)
        result = result + (a[i] * b[i]);

运行时库函数

线程环境相关

函数	作用
`omp_set_num_threads()`	设置使用的线程数
`omp_get_num_threads()`	获取使用的线程数
`omp_get_max_threads()`	获取可用最大的线程数
`omp_get_thread_num()`	获取线程编号
`omp_get_thread_limit()`	获取最大可用线程
`omp_get_num_procs()`	获取处理器数目

锁

函数	作用
`omp_init_lock()`	初始化锁
`omp_destroy_lock()`	删除锁
`omp_set_lock()`	加锁
`omp_unset_lock()`	解锁
`omp_test_lock()`	测试加锁
`omp_init_nest_lock()`	初始化嵌套锁
`omp_destroy_nest_lock()`	删除嵌套锁
`omp_set_nest_lock()`	加嵌套锁
`omp_unset_nest_lock()`	解嵌套锁
`omp_test_nest_lock()`	测试加嵌套锁

说明：

OpenMP 的锁是互斥锁，即 mutex

与 critical 指令的区别是， critical 指令保证同时只有一个线程运行内部的指令，而不是加锁

例子:

    #include <omp.h>

    omp_lock_t lock;
    omp_init_lock(&lock);
    int i = 0;
    #pragma omp parallel
    {
        while(true){
            omp_set_lock(&lock);
            x++;
            omp_unset_lock(&lock);
            break;
        }
    }
    omp_destroy_lock(&lock);

环境变量

变量	作用
`OMP_SCHEDULE`	`for` 指令的默认从句
`OMP_NUM_THREADS`	使用的最大线程数
`OMP_NESTED`	是否使用嵌套并行
`OMP_STACKSIZE`	设置子线程的栈大小
`OMP_WAIT_POLICY`	设置等待线程的行为
`OMP_MAX_ACTIVE_LEVELS`	设置嵌套并行的层级
`OMP_THREAD_LIMIT`	设置整个程序的最大线程数

SageMath

chimez

2021-02-19 14:27

量子力学算符代数计算

参见 Noncommutative Algebras in Sage.

例如要做满足 $[a, a^{\dagger}]=1, [b,b^{\dagger}]=1, [a,b]=0$ 的两个算符 $a,b$ 的代数计算

  F.<a, a_d, b, b_d> = FreeAlgebra(QQ)

  U_relations = {
      b_d * b: b * b_d - 1,
      a_d * a: a * a_d - 1,
  }

  U = F.g_algebra(U_relations)
  U._latex_names = ["a", r"a^\dagger", "b", r"b^\dagger"]
  U.inject_variables()

注意：

定义规则时，必须是后定义的字母乘在前面，例如定义时 F.<a,b,c> 规则中应该写作 c*b: ..., b*a: ... 等
变量的 latex 显示名顺序要和变量定义顺序相同
没有指定规则的乘法默认就是可对易的

安装为 jupyter 内核

jupyter kernelspec install ./SageMath/local/share/jupyter/kernels/sagemath

基本特性

直接输入数学运算公式可以作为计算器来用
继承自python的面向对象系统，即任何东西都是对象
在命令行模式下使用TAB补全，在notebook里似乎有bug
内置交互式帮助系统，在命令后加一个问号即可显示帮助 c.diameter?

符号运算

基础操作

`x = var('x')`	定义符号
`f = 1 - sin(x)^2`	定义一个函数
`print(f)`	普通打印函数
`print(maxima(f))`	美化打印函数
`f.simplify_trig()`	化简表达式
`f(x=pi/2)`	计算函数值
`ingtegrate(f, x)`	不定积分
`f.differentiate(2).substitude({x: 3/pi})`	计算二阶微分，并代替值
`a = var('a')`	x是默认隐含符号，需要更多符号就要定义
`S = solve(x^2 + x == a, x)`	解代数方程
`S[0].rhs()`	显示代数方程解
`show(plot(sin(x)+ sin(1.6*x)), 0, 40)`	解析画图
`F = factor(x^99 + y^99)`	分解因式
`F.expand()`	乘开因式

代数

用对象表示代数结构，群环域

组合数学

数值运算

基础操作

`3 + 5`	加法
`57.1 ^ 100`	乘方
`matrix([[1,2], [3,4]])`	创建矩阵，按行写
`g = sin(x) + (1- x^2)`	定义一个函数
`find_root(g, 0, 2)`	在区间[0, 2]寻根
`var('x y z')`	定义一些符号
`f = (1+ (y+ x^2)^2 + (1+ z+ y^2)^2)^2`	新定义一个函数
`minimize(f, [1, 2, 3], disp=1, algorithm='powell')`	计算最优化，指定算法和初值
`m = random_matrix(RDF, 500)`	创建复数500×500随机矩阵
`e = m.eigenvalues()`	计算矩阵的本征值
`w = [(i, abs(e[i])) for i in range(len(e))]`	循环方式创建list
`show(points(w))`	成对点画图
`factorial(100)`	计算阶乘
`N(pi, digits=100)`	保留100位小数，四舍五入
`z = Partitions(10^8).cardinality()`	计算数的划分并给出基数

包括以下几种浮点数类型：

Python	float, complex, decimal
SageMath specific	RDF, CDF, RQDF, CC, RR, RIF, CIF
included Systems	pari, maxima

绘图

可绘制的类型：

interact	动态交互图，可以显示一个滑块，滑动改变图
Regions Plot/ Contours	显示一个不等式区间
Density Plot	密度图，热力图
filled plot	可以上色显示两条线之间的区间
multiedge graph	可以显示一个有向图

这里有各种动态图的例子

技巧和常用代码片段

与 tex 联合使用

sagetex：可以在tex中显示sage命令的运算结果

首先需要安装 sagetex

注意：在网站下载的二进制包千万别手贱用 sage -i 自己再编译，由于包含文件不全，这样只会破坏已编译好的文件

    kpsewhich -var-value=TEXMFHOME                  # 获得TEXMFHOME目录位置，一般是$HOME/texmf
    cp -R SAGE_ROOT/local/share/texmf/tex TEXMFHOME # 把sagetex复制到目录里
    texhash                                         # 让latex重新索引包

使用方法

    \usepackage{sagetex}              % 引入包
    \usepackage[imagemagick]{sagetex} % 可选项
    \setlength{\sagetexindent}{10ex}  % 设置间隔
    $2+2=\sage{2+2}$                  % 行内计算sage表达式的数学结果
    \begin{sageblock}                 % 引用sage/python代码
      1+1
      var('a,b,c')
      eqn ~ [a+b*c==1, b-a*c==0, a+b=~5]
      s = solve(eqn, a,b,c)
    \end{sageblock} 
    \begin{sagesilent}      % 多行计算sage表达式
      e = 2
      e = 3*e + 1
    \end{sagesilent}
    $e=\sage{e}$            % 在同一文件内可以直接使用定义的变量
    \sageplot{E.plot(-3,3)} % sage画图
    \sagetexpause           % 在pause宏之间的sage表达式不被求值，便于修改
    \sagetexunpause

三角函数和指数之间转换

参考这个回答

   def exponentialize(x):
       """从三角函数形式变换成e指数形式"""
       return sageobj(x._maxima_().expand().exponentialize())

   def demoivre(x):
       """从e指数形式变换成三角函数形式"""
       return sageobj(x._maxima_().expand().demoivre())

C 指针与内存管理

chimez

2021-02-12 08:37

动态存储分配

内存分配函数

在 <stdlib.h> 头文件中的

malloc 分配内存，不初始化
calloc 分配内存，并清零
realloc 调制分配内存大小

内存分配函数返回 void * 类型，本质上只是内存地址

空指针

内存分配函数有可能返回空指针 NULL ，使用空指针会导致程序崩溃，所以要判断指针是否为空指针。

c语言中空指针为假，非空指针都为真。

常用如下方式进行判断

if (p == NULL) {
                ...;
}

malloc

原型 void *malloc(size_t size);

惯用法

double *p = malloc(N * sizeof(*p));

calloc

原型 void *calloc(size_t nmemb, size_t size); 分配 nmemb 个元素，每个元素大小为 size.

realloc

原型 void *realloc(void *ptr, size_t size);

扩展时，不会对新加进来的做初始化
扩展失败时，返回空指针并不影响原来的元素
传入空指针时，与 malloc 行为一致
第二个参数是0时，释放内存

一旦 realloc 返回，记得更新之前的所有指针，因为可能会把元素都复制到了别的地址

释放内存 free

原型 void free(void *ptr);

声明

声明说明符

存储类型。包括 auto, static, extern, register ，只能有一个，必须放在最前面
类型限定符。包括 const, volatile, restrict(C99) ，可以有一个或多个
类型说明符。包括基本类型、结构、枚举和联合
函数说明符。包括 inline(C99) ，只用于函数声明

解释复杂的声明符

始终从内往外读声明符
[]和()优先于 * ,即 *p[] 是数组， *f() 是函数

关键字

关键字	存储类别	存储期	作用域	链接	声明方式
可选 (`auto`)	自动	自动	块	无	块内
`register`	寄存器	自动	块	无	块内
无	静态外部链接	静态	文件	外部	函数外
`static`	静态内部链接	静态	文件	内部	函数外
`static`	静态无链接	静态	块	无	块内
`extern`	外部
`_Thread_local`
`const`	不可变
`inline`
`volatile`	代理，可由其它程序改变
`restrict`	限定唯一的指针
`_Atomic`

好的设计不应该使用文件作用域的变量
auto 关键字用于明确显示使用了与外部同名的变量名，与 c++ 中的完全不同，不建议使用
函数的默认类别是 extern 除非使用了 static 指定为模块私有
static 静态指的是生命周期一直在， const 才是不可变
const 变量只能初始化，不能修改
const float * pf 指针 pf 指向 const float 类型的对象，即对象内容不可变，但指针可以指向其它地址
float const * pf 同上
float * const pf 指针 pf 指向 float 类型的对象，即对象内容可变，但指针不能指向其它地址
即 const 在 * 前面表示指向的值不能变，在 * 后面表示指针本身不能变
函数的形参中有 const 表示传递的值不变
restrict 限制指针是指向内存的唯一和初始方式，用于优化

TODO

-Ofast -flto -march=native -funroll-loops -fPIC

pandas

chimez

2021-01-30 13:27

数据表格和矩阵之间的转换

从点的坐标表格 `(x,y,z)` 得到矩阵

   data = pd.DataFrame({"x": [1,1,2,2], "y":[3,4,3,4],"z":[9,8,7,6]})
   data_2x2 = data.pivot(index="x",columns="y", values="z")

   matrix = data_2x2.to_numpy()

从矩阵表格得到坐标表格

   data = pd.DataFrame({"x": [1,1,2,2], "y":[3,4,3,4],"z":[9,8,7,6]})
   data_2x2 = data.pivot(index="x",columns="y", values="z")

   d = data_2x2.unstack().reset_index()
   d.columns = ['y','x','z']

groupby: 分开-计算-合并

groupby guide 通过 groupby 可以做涉及到如下三个过程的一些操作：

根据某些标准将数据分类
对每类数据分别应用某个函数
将应用的结果合并回数据表中

分开数据后通常要应用三种函数

累积：计算每组的和、平均值、总数等
变换：进行与组内数据相关的计算并返回类似索引的列表，如对组内数据标准化、根据组的数据填上新的数据等
筛选：去掉某些不符合条件的组等

变换(transform)

   transformed = ts.groupby(lambda x: x.year).transform(
       lambda x: (x - x.mean()) / x.std()
   )

vimrc: vim

chimez

2021-01-28 18:33

正则表达式中使用非贪心匹配

.* 变成 .\{-} .+ 变成 .\{-1,\}

ffmpeg

chimez

2021-01-28 15:57

合并视频

  ffmpeg -f concat -i filelist.txt -c copy out.mkv

filelist.txt 内容如下，需要加 file ，还要用单引号括起来，或者加选项 -safe 0

file '/home/a1.mkv'
file '/home/a2.mkv'
...

转换格式

  ffmpeg -i in.mpeg -q:a 0 -q:v 0 out.mkv

选项	说明
`-q:a 0`	最高等级转换音频
`-q:v 0`	最高等级转换视频

分割视频

  ffmpeg -i in.mp4 -ss 00:10:00 -t 00:30:40 -acodec copy -vcodec copy output.mp4

命令的结果是分割得到00:10:00 –> 00:40:40 的原格式视频片段

选项	说明
`-ss 00:10:00`	从00:10:00开始分割
`-t 00:30:40`	分割时长00:30:40
`-acodec copy`	音频与原格式相同
`-vcodec copy`	视频与原格式相同

emacs

chimez

2021-01-26 18:28

从源码编译安装

查看特性

查看 system-configuration-features 中的信息

查看编译选项

查看 system-configuration-options 中的信息

   ./configure --without-compress-install --with-modules 'CFLAGS=-O2 -g3'

在正则表达式中使用lisp代码

使用 \, 后加一个表达式

  :s/\(.*?\)a\(.*?\)b/\1\,(+ (string-to-number \2) 100)/

正则表达式的非贪心匹配

non greedy regex 星号后加问号

:s/.*?//

文件编码问题

改变显示文件编码,解决乱码问题

使用命令

   revert-buffer-with-file-coding-system

修改编码

在按下 M-x 后使用命令改变编码再保存即可

   set-buffer-file-coding-system

常见的编码

中文
1. cp936
2. gb2312
3. gb18030
日文
1. cp932

对齐排版代码

使用 align-regexp, sort-regexp-fields, sort-columns

Tramp 远程编辑

配置 ssh

注意 tramp 的远程登陆只支持 bash, 不支持 zsh, 所以要进行设置，tramp 登陆用 namebash, ssh 登陆用 namezsh 就行了

文件 $HOME/.ssh/config

   Host namebash
       HostName 10.0.0.1
       User username
       RequestTTY yes
       RemoteCommand bash

   Host namezsh
       HostName 10.0.0.1
       User username
       RequestTTY yes
       RemoteCommand zsh

Windows 相关问题

创建快捷方式

实现以下功能

隐藏终端窗口
启用 conda 环境
设置默认工作路径

   $WshShell = New-Object -comObject WScript.Shell
   $Shortcut = $WshShell.CreateShortcut("$([environment]::GetFolderPath("Desktop"))\emacs.lnk")
   $Shortcut.TargetPath = "PowerShell.exe"
   $Shortcut.Arguments = " -Windowstyle Hidden -Command ""& 'C:\miniconda3\shell\condabin\conda-hook.ps1' ; conda activate 'C:\miniconda3'; conda activate py39; C:\msys64\mingw64\bin\runemacs.exe"""
   $Shortcut.WorkingDirectory = "C:\"
   $Shortcut.Save()

vimrc: vim配置文件

chimez

2021-01-24 13:57

<src

Write your post here.

PyTables: 分层数据集

安装

基本使用

创建一个 HDF5 文件

在文件中创建任意层的数据集

在数据集中创建若干个数据表

遇到的问题

在不支持 flock 的集群上运行

h5py: 储存数据

存

取

OpenMP 共享内存并行

参考资料

编译器指令

parallel 指令

for 指令

sections, section 指令

single 指令

task 指令

master 指令

critical 指令

barrier 指令

taskwait 指令

atomic 指令

flush 指令

ordered 指令

threadprivate 指令

reduction 从句

运行时库函数

线程环境相关

锁

环境变量

量子力学算符代数计算

安装为 jupyter 内核

基本特性

符号运算

基础操作

代数

组合数学

数值运算

基础操作

绘图

技巧和常用代码片段

与 tex 联合使用

sagetex：可以在tex中显示sage命令的运算结果

三角函数和指数之间转换

动态存储分配

内存分配函数

空指针

malloc

calloc

realloc

释放内存 free

声明

声明说明符

解释复杂的声明符

关键字

TODO

数据表格和矩阵之间的转换

从点的坐标表格 (x,y,z) 得到矩阵

从矩阵表格得到坐标表格

groupby: 分开-计算-合并

变换(transform)

正则表达式中使用非贪心匹配

合并视频

转换格式

分割视频

从源码编译安装

查看特性

查看编译选项

在正则表达式中使用lisp代码

正则表达式的非贪心匹配

文件编码问题

改变显示文件编码,解决乱码问题

修改编码

常见的编码

对齐排版代码

Tramp 远程编辑

配置 ssh

Windows 相关问题

在不支持 `flock` 的集群上运行

从点的坐标表格 `(x,y,z)` 得到矩阵