Skip to content

Welcome!

My flagship website has moved to the new address. This site is only for documentation purpose.

Memo of several opt-used Python commands

  • Installation
  • To install python packages to specified directory, two options are available:

    easy_install --install-dir=/your/directory package
    pip install --install-option="--install-purelib=/your/directory" package

    Uninstall python packages:
    pip uninstall package

    Search python package (for example, scikit-learn):
    pip search scikit

    Generate a profile of python code by kernprof, then analyze the profile by pstats:
    python -m kernprof pyfile
    python -m pstats pyfile.prof

  • Useful commands in IPython
  • Reload a module (minres, for example):
    import minres as m
    reload(m.MINRES)

    Display all variables in working space
    who

    Clear variables in current working space
    reset -f

    Run Linux command from Python
    os.system(cmd)

    Return the position of the largest element of a list
    L.index(max(L))

    Check the class to which a variable belongs
    type(var)

    Execute a script file:
    execfile(filename)

  • Editing
  • Hide Show mode in Emacs:

    C-c @ C-M-s show all
    C-c @ C-M-h hide all
    C-c @ C-s show block
    C-c @ C-h hide block
    C-c @ C-c toggle hide/show

    Some inequalities on block matrix inverse norm

    Suppose \(H\in\mathbb R^{n\times n}\) is a positive-definite matrix. Denote the index set \(\mathcal I, J\subseteq \{1,2,\ldots,n\}\), and \(H_{\mathcal I\mathcal J}\) as the submatrix of \(H\) containing the rows of \(\mathcal I\) and columns of \( \mathcal J\). An interesting question arise as whether \[ \|(H_{\mathcal I\mathcal I})^{-1}\|\leq \| (H^{-1})_{\mathcal I\mathcal I}\|\leq \| H^{-1}\|\] for arbitrary matrix norm?

    By applying matrix inversion lemma on blockwise inversion matrix, \(\|(H_{\mathcal I\mathcal I})^{-1}\|_2\leq\|(H^{-1})_{\mathcal I\mathcal I}\|_2\). Therefore, it is easily seen that \[ \|(H_{\mathcal I\mathcal I})^{-1}\|_2\leq \| (H^{-1})_{\mathcal I\mathcal I}\|_2\leq \| H^{-1}\|_2\mbox{ is correct.}\]

    However, both inequalities may fall apart for other matrix norms, say \(\|\cdot\|_1\) or \(\|\cdot\|_\infty\).
    \[
    H =
    \begin{pmatrix}
    \frac{7}{6} & \frac{1}{6} & \frac{29}{12}\\
    \frac{1}{6} & \frac{1}{6} & \frac{5}{12}\\
    \frac{29}{12} & \frac{5}{12} & \frac{143}{24}
    \end{pmatrix},
    \quad \mathcal I = \{1,2\}.
    \]

    It is verifiable that for this example \( \| (H^{-1})_{\mathcal I\mathcal I}\|_\infty \leq\| H^{-1}\|_\infty < \|(H_{\mathcal I\mathcal I})^{-1}\|_\infty\).

    More intriguingly, suppose \(H\) is an \(M\)-matrix, i.e. positive definite with only nonpositive off-diagonal entries. It can be shown that \( 0\leq(H_{\mathcal I\mathcal I})^{-1}\leq (H^{-1})_{\mathcal I\mathcal I} \), and therefore \(\|(H_{\mathcal I\mathcal I})^{-1}\|_\infty\leq \| (H^{-1})_{\mathcal I\mathcal I}\|_\infty \leq \| H^{-1}\|_\infty\).

    Counterexamples in plain primal-dual active-set method

    Consider the simple bound-constrained QP of the form
    \[\min_{x\in\mathbb R^n}\{\frac{1}{2}x^THx + c^Tx:\; x\leq u\}.\]
    Suppose \(H\succ 0\) which guarantees that the primal-dual active-set algorithm is well-defined

    1. Input an initial partition \( (\mathcal A,\mathcal I) \)
    2. Loop
    3. \(\quad\) Set \(x_{\mathcal A}\gets u_{\mathcal A},\ z_{\mathcal I}\gets 0\), compute \(x_{\mathcal I}\) by solving
      \[H_{\mathcal I\mathcal I}x_{\mathcal I} = -H_{\mathcal I\mathcal A}u_{\mathcal A} -c_{\mathcal I},\]
      \(\quad\;\)and set \( z_{\mathcal A} = -H_{\mathcal A:}x -c_{\mathcal A}\).
    4. \(\quad\) Let \(V_P:=\{i\in\mathcal I:\ x_i > u_i\}\) and \(V_D:=\{i\in\mathcal A:\ z_i < 0\}\).
    5. \(\quad\) If \(V_P\cup V_D=\emptyset\), then break and return \( (x,z) \).
    6. \(\quad\) Set \(\mathcal A\gets (\mathcal A\backslash V_D)\cup\mathcal V_P\), and \(\mathcal I\gets (\mathcal I\backslash V_P)\cup\mathcal V_D\).
  • An example that fails to converge
  • PDAS is a very simple and efficient algorithm when \(H\) is an \(M\)-matrix, i.e. positive definite and with nonpositive off-diagonal entries. However, simply positive definiteness may result to cycles for \(n\geq 3\). Consider the following example
    \[
    H=\begin{pmatrix}
    4 & 5 & -5\\
    5 & 9 & -5\\
    -5&-5 & 7
    \end{pmatrix}, c=\begin{pmatrix}
    2\\ 1\\ -3
    \end{pmatrix}, u=\begin{pmatrix}
    0\\ 0 \\ 0
    \end{pmatrix}.
    \]
    For convenience, employ a tuple \( (i_1,i_2,i_3) \) to indicate to which set a corresponding index belongs. For example \( (i,i,u) \) means 1,2 belongs to inactive-set \(\mathcal I\) and 3 belongs to active-set \(\mathcal A\). For the specified example, the optimal partition is with tuple \(i,u,u\). The PDAS algorithm is essentially iterating from one partition to another and the following figure shows the iteration when PDAS is applied.

    Nonconvergence of PDAS

    The figure shows that 6 of the possible 8 partitions will cycle and can never reach the optimal one!

  • An example that impedes inexact linear equations solves
  • Step 3 of PDAS requires solving linear equations to obtain \(x_{\mathcal I}\). To save computation, it is natural to consider solving the linear equations inexactly, generate an inexact \(\tilde x_{\mathcal I}\) with residual \(r\). As for \(V_P\) and \(V_D\), define \(\tilde V_P\) and \(\tilde V_D\) from \(\tilde x\) and \(\tilde z\). However, it may happen that requiring \(\tilde V_P\subseteq V_P\) and \(\tilde V_D\subseteq V_D\) may enforce \(r=\boldsymbol 0\). In other words, we can not identify subsets of the violations \(V_P\) and \(V_D\) unless the exact solution \((x,z)\) is known. The counterexample is
    \[H=\begin{pmatrix}
    2 & 0 & 0\\
    0 & 2 &-1\\
    0 &-1 & 2
    \end{pmatrix}, c =\begin{pmatrix}
    1\\ -1 \\ -1
    \end{pmatrix}, u = \boldsymbol 1, \mathcal A=\{1,2\}, \mathcal I=\{3\}.
    \]
     

    An application of Hölder’s inequality


    One way of describing Hölder’s inequality is this: Let \(x,y\in\mathbb R^n\), \(\|\cdot\|_p\) and \(\|\cdot\|_q\) are pre-norm and dual-norm respectively, then \(|x^Ty|\leq \|x\|_p \|y\|_q\). Note taht given a pre-norm \(f(x)\), the dual norm is defined as \( f^D(y):=\max_{f(x)=1} y^*x \), in our case \(\frac{1}{p} + \frac{1}{q} = 1\). Below is one application that is vital in one of our algorithm designs.

    Proposition: Let \(A\in\mathbb R^{m\times n}\) and \(x\in\mathbb R^n\), \(\|\cdot\|_p\) and \(\|\cdot\|_q\) are pre-norm and dual-norm respectively, then

    \[|Ax|\leq \| A^T\|_p \| x\|_q\boldsymbol 1.\]

    Proof: Denote \(A^T = (a_1,a_2,\ldots,a_m)\) where \(a_i^T\) is the \(i\)th row of \(A\). Then by Hölder’s inequality, we have \(|a_i^Tx|\leq \|a_i\|_p\|x\|_q\). To be concrete, the matrix-vector multiplication

    \[ |Ax| =\begin{pmatrix}
    |a_1^Tx|\\
    |a_2^Tx|\\
    \vdots\\
    |a_m^Tx|
    \end{pmatrix}
    \leq \begin{pmatrix}
    \|a_1\|_p\|x\|_q\\
    \|a_2\|_p\|x\|_q\\
    \vdots\\
    \|a_m\|_p\|x\|_q
    \end{pmatrix}
    = \begin{pmatrix}
    \|a_1\|_p\\
    \|a_2\|_p\\
    \vdots\\
    \|a_m\|_p
    \end{pmatrix} \|x\|_q.
    \]
    Next, the fact that \(a_i^T\) is a row of \(A\) thus \(a_i\) a column of \(A^T\) implies that \(\|a_i\|_p= \|A^T\boldsymbol e_i\|_p \leq \|A^T\|_p\). Therefore, we further have
    \[
    |Ax|\leq \|A^T\|_p\|x\|_q\boldsymbol 1.
    \]

    Remark: In our algorithm, we need to bound a matrix-vector multiplication term \(H^{-1}r\) where we do not want to compute \(H^{-1}\). The proposition shows that \(|H^{-1}r|\leq \|H^{-T}\|_p \|r\|_q\boldsymbol 1\) where an upper bound of \(\|H^{-T}\|_p\) sometimes is cheaper to obtain.

    Expression of emotion: Lehigh News in the decade 2002 – 2013


    This article is inspired by The Expression of Emotions in 20th Century Books, which crunches huge amount of data from Google books to illustrate the effectiveness of simple statistics by leveraging big data, and AFINN, a list of words  with mood values betwee -5 (extremely negative) and 5 (extremely positive).  The objective is to analyze the mood changes of the decade in Lehigh University. The data is obtained by crawling Lehigh News Archive links, which is doable from Python.

    • Results

    The following figure shows the annual emotional fluctuation between 2002 and 2013. The red curve represents the intensity of positive emotions (positive value), the blue curve of negative emotions (negative value), and the black curve the average to show the overall mood score. The bar chart reflects the absolute value of positive and negative moods.

    The figure shows that the year 2004, with the lowest mood score, is gloomy. Looking back at the news posted in 2004, there are such news as In memoriam: Muriel Mallin Berman, and Mountain Hawks dream season ends in the first half of the year. The mood score increased almost steadily from 2005 to 2008, then dipped in 2009: the subprime mortgage crisis? Probably. The years of 2010 through 2011 looks no better than 2009. On November 2011, many negative news are posted on Health Advisory. In 2012, Lehigh has undergone sandy A “perfect storm” rages through campus, but also experienced the joy of victory over Duke in NCAA Tourney.

    Lehigh Mood in the Decade

     

    • Methodology

    The AFINN list includes 2477 words with values between -5 and -5. Denote the value of the \(i\)th word as \(v_i\), we count the frequency \( (c_1,c_2,\ldots,c_n) \) for each of the 2477 words for each year. The overall mood score of each year is given by
    \[\sum_{i=1}^n \frac{c_iv_i}{c_{the}},\]
    where \(c_{the}\) is the frequency of the word “the”, widely used as a measure of the length of an article. The intensity of positive/negative mood is computed by summing on those \(i\)s with positive/negative \(v_i\).

    • Remark

    It may be more interesting to show the monthly mood value fluctuations which is possibly a better indicator of the mood. However, this will entail the web crawler to collect news of each month, which is hard because of the way Lehigh News Archive is organized.

    Algorithm Series 2: Longest Substring Without Repeating Characters

    This problem derives from the “Longest Substring Without Adjacent Repeating Characters“, the difference lies on how “longest” is defined.

    Given a string, find the length of the longest substring without repeating characters. For example, the value for “abcabcbb” is “abc”, which the length is 3. For “bbbbb” the longest substring is “b”, with the length of 1.

    The algorithm looks slightly simpler:

    def lengthOfLongestSubstring2(S):
        '''
        Given a string, find the length of the longest substring without repeating characters.
        For example, the longest substring without repeating letters for "abcabcbb" is "abc" with length 3.
        For "bbbbb" the longest substring is "b", with the length of 1.
        '''
        # Base cases
        lenS = len(S)
        if lenS <= 1:
            return lenS
        # Recursively obtain lengthOfLongestSubstring2(S-1)
        lenCand = lengthOfLongestSubstring2(S[:-1])
        lenCand = max(lenCand, len(set(S[-lenCand-1:])) )
        return lenCand
    

    And again, the complexity is O(n^2) by previous analysis.

    Algorithm Series 1: Longest Substring Without Adjacent Repeating Characters

    The elegance of Python decides that it is ideal for demonstrating algorithms. This example illustrates a recursive algorithm to solve the “Longest Substring without Adjacent Repeating Characters” problem:

    Given a string, find the length of the longest substring without adjacent repeating characters. For example, the value for “abcabcbb” is “abcabcb”, which the length is 7. For “bbbbb” the longest substring is “b”, with the length of 1.

    The Python code is quite self-explanatory:

    def lengthOfLongestSubstring(S):
        '''
        Given a string, find the length of the longest substring without repeating characters.
        For example, the value for "abcabcbb" is "abc", which the length is 3.
        For "bbbbb" the longest substring is "b", with the length of 1.
        '''
        # Base cases
        lenS = len(S)
        if lenS <= 1:
            return lenS
        # Recursively obtain lengthOfLongestSubstring(S-1)
        lenCand = lengthOfLongestSubstring(S[:-1])
    
        if S[-1] == S[-2]:
            return lenCand
        elif is_full(S[-lenCand-1:-1]):
            return lenCand + 1
        else:
            return lenCand
    
    def is_full(S):
        '''
        A function to help decide if the string S itself is the longest substring without repeating characters.
        '''
        if len(S) <= 1:
            return True
        for i in range(len(S)-1):
            if S[i] == S[i+1]:
                return False
    
        return True
    
    if __name__ == "__main__":
    
        assert lengthOfLongestSubstring([]) == 0, 'Empty list should return 0.'
        assert lengthOfLongestSubstring('a') == 1, 'One character should return 1.'
        assert lengthOfLongestSubstring('a'*6) == 1, 'Should return 1.'
        assert lengthOfLongestSubstring('abbababaabb') == 6, 'Should return 6.'
        S = "wteubifcbeajzhnzvdrxyismtdgbscxqyclzksdnwgzypmxlsqisaceuglvapurnyepkwuavaztqnsbhjlzjoefurcw"
        assert lengthOfLongestSubstring(S) == 91, 'A random string that returns 91.'

    In short, to obtain lengthOfLongestSubstring(S), we need to know the value of the substring of S, where the substring is S with the tail element removed. Therefore, we have O(n) = O(n-1) + |candLen| and the algorithm is O(n^2). Dynamic programming may yield a more efficient algorithm, but the word “recursive” come up from my mind at first glance of this problem.

    A quick tutorial for installing CUTEr, an optimization test set

    According to the official introduction, CUTEr is the acronym for:

    The Constrained and Unconstrained Testing Environment, revisited/safe threads, (CUTErCUTEst) contains a set of test problems for linear and nonlinear optimization.

    This article serves as a quick-start tutorial for installing CUTEr. Crucial steps are:

    1. Download necessary packages.
    2. Install cuter2 and sifdec2.
    3. Setup environmental variables.
    4. Setup Matlab.
    5. Test installation.

    Throughout this article, assume the installation is under CUTErInstallation/ (absolute path: $HOME/CUTErInstallation/) which contains two subdirectories cuter2/ and sifdec2/ .

    1.  Download packages.

    • Change working directory to CUTErInstallation
    $ cd CUTErInstallation
    • Download relevant files to the corresponding folders
    $ svn co http://tracsvn.mathappl.polymtl.ca/SVN/cuter/sifdec/branches/SifDec2 ./sifdec2
    $ svn co http://tracsvn.mathappl.polymtl.ca/SVN/cuter/cuter/branches/cuter64 ./cuter2

    or for 32 bit machine:

    $ svn co http://tracsvn.mathappl.polymtl.ca/SVN/cuter/cuter/branches/CUTEr2 ./cuter2
    • Download and unzip mastsif problems
    $ wget ftp://ftp.numerical.rl.ac.uk/pub/cute/mastsif_small.tar.gz
    $ wget ftp://ftp.numerical.rl.ac.uk/pub/cuter/mastsif_large.tar.gz
    $ tar -zvxf mastsif_large.tar.gz
    $ tar -zvxf small.tar.gz

    2. Install cuter2 and sifdec2.

    • Install cuter
    $ cd cuter2/
    $ ./install_cuter

    Select Platform([5] PC), Fortran compiler ([7] gfortran), C compiler ([3] GNU g++), install precision ([D] Double), install size ([L] Large) and Y for all other options.

    • Install sifdec
    $ cd ../sifdec2/
    $ ./install_sifdec

    Select Platform([5] PC), Fortran compiler ([7] gfortran), install precision ([D] Double), install size ([L] Large) and Y for all other options.

    3.  Set up environmental variables (depending on step 2 of your selection of settings).

    • Open the initialization file
    $ Emacs $HOME/.bashrc
    • Add the following lines:
    export SIFDEC="$HOME/CUTErInstallation/sifdec2"
    export MYSIFDEC="$HOME/CUTErInstallation/sifdec2/SifDec.large.pc.lnx.gfo"
    
    export MANPATH="$MANPATH:$HOME/CUTErInstallation/sifdec2/common/man"
    export PATH="$PATH:$HOME/CUTErInstallation/sifdec2/SifDec.large.pc.lnx.gfo/bin"
    
    export CUTER="$HOME/CUTErInstallation/cuter2"
    export MYCUTER="$HOME/CUTErInstallationcuter2/CUTEr.large.pc.lnx.gfo"
    
    export MANPATH="$MANPATH:$HOME/CUTErInstallation/cuter2/common/man"
    export PATH="$PATH:$HOME/CUTErInstallation/cuter2/CUTEr.large.pc.lnx.gfo/bin"
    export MATLABPATH="$MATLABPATH:$HOME/CUTErInstallation/cuter2/CUTEr.large.pc.lnx.gfo/double/bin:$HOME/CUTErInstallation//cuter2/common/src/matlab"
    export LD_LIBRARY_PATH="$LD_LIBRARY_PATH:$HOME/CUTErInstallation/cuter2/CUTEr.large.pc.lnx.gfo/double/lib"
    
    export MASTSIF="$HOME/CUTErInstallation/mastsif"

    Reload the bash file:

    $ source $HOME/.bashrc

    4. Setup Matlab.

    • Add the following directory to Matlab working path
    $HOME/CUTErInstallation/cuter2/CUTEr.large.pc.lnx.gfo/double/bin/
    $HOME/CUTErInstallation/cuter2/common/scr/matlab/

    and setup MEX in Matlab:

    >>mex -setup

    and choose meshopts.sh.

    5. Now cuter problems should be ready for test, commands

    $ runcuter --package mx --decode GENHS28

    , if CUTEr is installed successfully, should prompt:

    sifdecode  GENHS28
    
     Problem name: GENHS28   
    
     Double precision version will be formed. 
    
     The objective function uses        2 linear groups
    
     There are        8 linear equality constraints
    
     There are       10 free variables

    and in Matlab:

     >> prob = cuter_setup()

    should yield:

    prob = 
    
             n: 10
             m: 8
          nnzh: 19
          nnzj: 24
             x: [10x1 double]
            bl: [10x1 double]
            bu: [10x1 double]
             v: [8x1 double]
            cl: [8x1 double]
            cu: [8x1 double]
        equatn: [8x1 logical]
        linear: [8x1 logical]
          name: 'GENHS28   '

    Reference: http://tracsvn.mathappl.polymtl.ca/trac/cuter

    假期回国科研必备利器

    • Ubuntu/Linux Mint软件源

    以LM14(Nadia) 为例,默认的软件源(/etc/apt/sources.list)是

    deb http://packages.linuxmint.com/ nadia main upstream import
    deb http://archive.ubuntu.com/ubuntu/ quantal main restricted universe multiverse
    deb http://archive.ubuntu.com/ubuntu/ quantal-updates main restricted universe multiverse
    deb http://security.ubuntu.com/ubuntu/ quantal-security main restricted universe multiverse
    deb http://archive.canonical.com/ubuntu/ quantal partner
    deb http://packages.medibuntu.org/ quantal free non-free

    # deb http://archive.getdeb.net/ubuntu quantal-getdeb apps
    # deb http://archive.getdeb.net/ubuntu quantal-getdeb games

    从默认软件源下载texlive和auctex可以用龟速来形容。可以看出,http://packages.linuxmint.com/ 是LM默认的软件源, http://archive.ubuntu.com/ubuntu/等是Ubuntu的默认软件源。替代以相对应的国内软件源可以极大提高软件下载速度:

    可以选择的LM的国内软件源

    http://mirrors.ustc.edu.cn/linuxmint/(中科大)

    http://jx.tju.zyrj.org/linuxmint/(天大)

    Ubuntu的软件源

    http://mirrors.163.com/ubuntu(速度还不错)

    http://tw.archive.ubuntu.com/ubuntu/

    http://mirror.bjtu.edu.cn/ubuntu/

    http://ftp.sjtu.edu.cn/ubuntu/

    http://mirror9.net9.org/ubuntu/

    http://ubuntu.csie.ntu.edu.tw/ubuntu/

    http://ubuntu.xjtuns.cn/ubuntu/

    • Lehigh的VPN配置

    Lehigh的VPN对访问学校资源很必要,无奈http://sslvpn.lehigh.edu/不稳定,浏览器安装cisco的AnyConnectVPN经常失败。Alternative方法是安装 openconnect 和 network-manager-openconnect

    sudo apt-get install openconnect network-manager-openconnect

    比AnyConnectVPN好用太多了。安装完即可在Network Setting中设置VPN,好处不必多说,访问Lehigh的Hdrive和系里的Server是肯定没问题了,上Gmail也不再掉线了。

    • Dropbox配置

    安装客户端的话需要按这里下载一个东西替代/home文件夹下的.dropbox-dist文件夹,更简单的做法是找个Dropbox的替代品,云诺貌似模仿的最好支持的操作系统最多,32位系统需要先安装ia32-libs。

    • Google加速

    需要下载一个list(google.go), 添加到host中:

    cat google.go >> /etc/hosts

    详见暗度陈仓:修改Hosts改善谷歌浏览速度

    AUCTex and RefTex in Emacs

    AUCTex has many awesome features in editing LaTex files within Emacs: preview, outline, folding display, among others. RefTex is a nice tool to manage cross references, in an elegant way. The price we pay for those conveniences is simply a configuration of the .emacs file. The following is my awkward example by revision from many other resources but it works well enough for me:

    (setq Tex-PDF-mode t)
    (load "auctex.el" nil t t)
    (require 'tex-mik)
    (load "preview-latex.el" nil t t)
    
    (setq Tex-save-query nil)
    (setq TeX-auto-save t)
    (setq TeX-parse-self t)
    (setq-default TeX-master nil)
    (setq Tex-electric-sub-and-superscript 1)
    (setq reftex-plug-into-AUCTeX t)
    
    (add-hook 'LaTeX-mode-hook 'visual-line-mode)
    (add-hook 'LaTeX-mode-hook 'flyspell-mode)
    (add-hook 'LaTeX-mode-hook 'LaTeX-math-mode)
    (add-hook 'LaTeX-mode-hook 'turn-on-reftex)
    (add-hook 'LaTeX-mode-hook
              (lambda ()
                (LaTeX-add-environments
                  '("theorem" LaTeX-env-label)
                  '("lemma" LaTeX-env-label)
                  '("proof" LaTeX-env-label))))
    
    (custom-set-variables
      ;; custom-set-variables was added by Custom.
      ;; If you edit it by hand, you could mess it up, so be careful.
      ;; Your init file should contain only one such instance.
      ;; If there is more than one, they won't work right.
     '(preview-gs-options (quote ("-q" "-dNOPAUSE" "-DNOPLATFONTS" "-dPrinted" "-dTextAlphaBits=4" "-dGraphicsAlphaBits=4"))))

    Outline:
    This feature allows you to view the big picture of the article by only displaying the sections, subsections etc.. To enable outline, turn on outline mode:

    M-x outline-minor-mode

    Next, C-c @ C-t outlines all sections and subsections:

    Move cursor to the section/subsection to be edited, C-c @ C-e displays the section/subsection:

    Additionally, C-c C-p C-b preview all math formulas in the buffer:

    Move closer to the citation, C-c ] guides to search for a desired reference:

    A reference card should help a lot for quick start.