物理の駅 Physics station by 現役研究者

テクノロジーは共有されてこそ栄える

Visual Studio: Cで書かれたライブラリをC++から呼ぶ方法

C言語で書かれたライブラリ (.libファイル) を C++からは直接呼ぶことができない。検証してみる。

libtestというプロジェクトで、構成の種類をスタティックライブラリ (.lib)にしたプロジェクトを作成した

C言語で書かれたライブラリと、サンプルプログラムを以下のように作成。

libtest.c

#include "libtest.h"

double libtest()
{
    return 0.0;
}

libtest.h

#pragma once

double libtest();

simple.cpp

#include <iostream>
#include <libtest.h>

int main() {
    std::cout << libtest() << std::endl;

    return 0;
}

これで、libtestとsimpleをビルドすると、simpleの方のビルドは失敗した。出力は以下の通り。

ビルドを開始しました...
1>------ ビルド開始: プロジェクト: libtest, 構成: Release Win32 ------
1>libtest.c
1>libtest.vcxproj -> C:\Users\Masahiro\source\repos\libtest\Release\libtest.lib
2>------ ビルド開始: プロジェクト: simple, 構成: Release Win32 ------
2>simple.obj : error LNK2001: 外部シンボル "double __cdecl libtest(void)" (?libtest@@YANXZ) は未解決です
2>C:\Users\Masahiro\source\repos\libtest\Release\simple.exe : fatal error LNK1120: 1 件の未解決の外部参照
2>プロジェクト "simple.vcxproj" のビルドが終了しました -- 失敗。
========== ビルド: 1 正常終了、1 失敗、0 更新不要、0 スキップ ==========

Visual C++ユーザーにはおなじみの 外部シンボル は未解決です である。libtest.cには確かに関数の定義が書かれているのに、呼ばれていない。

libtest.c のファイル名、libtest.cpp (拡張子 c→cpp)に変更して、libtestとsimpleプロジェクトをビルドする。

ビルドは成功した。出力は以下の通り。

ビルドを開始しました...
1>------ ビルド開始: プロジェクト: libtest, 構成: Release Win32 ------
1>libtest.cpp
1>libtest.vcxproj -> C:\Users\Masahiro\source\repos\libtest\Release\libtest.lib
2>------ ビルド開始: プロジェクト: simple, 構成: Release Win32 ------
2>コード生成しています。
2>0 of 3 functions ( 0.0%) were compiled, the rest were copied from previous compilation.
2>  0 functions were new in current compilation
2>  0 functions had inline decision re-evaluated but remain unchanged
2>コード生成が終了しました。
2>simple.vcxproj -> C:\Users\Masahiro\source\repos\libtest\Release\simple.exe
========== ビルド: 2 正常終了、0 失敗、0 更新不要、0 スキップ ==========

次に、libtest.cpp を libtest.cに戻して、simple.cpp のヘッダーの宣言を以下のように変えてみる。

extern "C" {
#include <libtest.h>
}

ビルドしたときの出力。

ビルドを開始しました...
1>------ ビルド開始: プロジェクト: libtest, 構成: Release Win32 ------
1>libtest.c
1>libtest.vcxproj -> C:\Users\Masahiro\source\repos\libtest\Release\libtest.lib
2>------ ビルド開始: プロジェクト: simple, 構成: Release Win32 ------
2>コード生成しています。
2>Previous IPDB and IOBJ mismatch, fall back to full compilation.
2>All 3 functions were compiled because no usable IPDB/IOBJ from previous compilation was found.
2>コード生成が終了しました。
2>simple.vcxproj -> C:\Users\Masahiro\source\repos\libtest\Release\simple.exe
========== ビルド: 2 正常終了、0 失敗、0 更新不要、0 スキップ ==========

成功した。

下記の参照リンクによると

結局のところ、C++からCのモジュールを呼び出す際にはextern "C"を用いるということです。なぜなら、C++コンパイラでは、マングリングによって、シンボル名が関数名ではなくなり、Cコンパイラとの互換性がなくなってしまうからです。

結論: C++から、C言語で書かれたライブラリを呼ぶとき、ヘッダー宣言時に、extern "C"{} で囲もう。

参照

qiita.com

Pythonから外部プログラムを実行して標準出力をstrで得る

WindowsでもC++を直接呼ぶ方法があるらしいが、めんどくさいので引数で情報を渡して、標準出力で結果を得ることにする。

import subprocess

exe = 'puroguramu.exe'
arg = "hikisu" #引数が1つの例

proc = subprocess.Popen([exe, arg], stdout=subprocess.PIPE)

for line in proc.stdout: #1行ずつ標準出力を得る
    print(line.decode('utf-8')) # byteからutf-8にデコード

proc.wait() #プロセスが終わるまで待つ
print(proc.returncode) #戻り値

参照

blog.imind.jp

Python lmfitを使って、マルチガウシアンでフィッティングする

phst.hateblo.jp

の続編

mzfit が使っている tensorflow が重すぎなのと、モジュール読み込みに時間がかかるのと、機能的に不足する点があったのと、マルチプラットフォームであるべきという思想に反するので、lmfitを使ってマルチガウシアンでのフィッティングを書いてみた。

lmfitの公式ドキュメントのFitting Multiple Peaks用コードを参考にした。

データ生成部分。

import numpy as np

np.random.seed(seed=0)
signal1 = np.random.normal(5, 1, 2000)
signal2 = np.random.normal(10, 2, 10000)
signal3 = np.random.normal(15, 0.5, 5000)
data = np.concatenate([signal1, signal2, signal3])

hist, bin_edges = np.histogram(data,bins=100,range=(2,18))
bin_center = (bin_edges[1:]+bin_edges[:-1])*0.5
arr_x = []
arr_y = []
arr_yerror = []
for x,y in zip(bin_center,hist):
    if y==0:continue
    arr_x.append(float(x))
    arr_y.append(float(y))
    arr_yerror.append(y**0.5)
arr_x = np.array(arr_x)
arr_y = np.array(arr_y)
arr_yerror = np.array(arr_yerror)

3つのガウシアンモデルで関数形を用意し、初期値を与える。初期値を全く与えなくてもfittingはうまくいったが、作法として初期値は与えておくべきである。set関数には、valueだけでなく、パラメータの範囲はminmaxとして引数を与えることができる。

from lmfit.models import ExponentialModel, GaussianModel

gauss1 = GaussianModel(prefix='g1_')
pars = gauss1.make_params()
gauss2 = GaussianModel(prefix='g2_')
pars.update(gauss2.make_params())
gauss3 = GaussianModel(prefix='g3_')
pars.update(gauss3.make_params())

pars['g1_center'].set(value=3)
pars['g1_amplitude'].set(value=100)
pars['g2_center'].set(value=10)
pars['g2_sigma'].set(value=1)
pars['g2_amplitude'].set(value=100)
pars['g3_center'].set(value=16)
pars['g3_sigma'].set(value=1)
pars['g3_amplitude'].set(value=100)

mod = gauss1 + gauss2 + gauss3

evalで初期パラメータを使ったarr_yの値を取得する。fitでフィッティングを行う。先の記事のとおり、分散共分散行列を正しく求めるため、引数にweights=1/arr_yerrorscale_covar=False を与える。

init = mod.eval(pars, x=arr_x)
result = mod.fit(arr_y, pars, x=arr_x, weights=1/arr_yerror, scale_covar=False)

result.result.params

グラフ化

import matplotlib.pyplot as plt

fig, axes = plt.subplots(1, 2, figsize=(12.8, 4.8))
axes[0].errorbar(arr_x, arr_y, yerr=arr_yerror, marker='.', linestyle="")
axes[0].plot(arr_x, init, '--', label='initial fit')
axes[0].plot(arr_x, result.best_fit, '-', label='best fit')
axes[0].legend()

comps = result.eval_components(x=arr_x)
axes[1].errorbar(arr_x, arr_y, yerr=arr_yerror, marker='.', linestyle="")
axes[1].plot(arr_x, comps['g1_'], '--', label='Gaussian component 1')
axes[1].plot(arr_x, comps['g2_'], '-', label='Gaussian component 2')
axes[1].plot(arr_x, comps['g3_'], ':', label='Gaussian component 3')
axes[1].legend()

plt.show()

マルチガウシアンでフィット

Python ウェブサイトで公開されているネットワーク障害の頻度を解析する

名古屋大学の障害情報一覧にあるネットワーク障害の発生年月日、発生月をPythonで解析してみる。2022/05/14 時点でのHTML構造に対応しているので、ホームページの構造が更新されると使えなくなることをご承知おきください。

データをダウンロード

import urllib.request
import datetime

suffix_urls = [""]
for i in range(2,14):
    suffix_urls.append(f"page{i}/")
datetimes = []
months = []
for suffix_url in suffix_urls:
    url = f"https://icts.nagoya-u.ac.jp/ja/information/trouble/{suffix_url}"
    filename = "a.html"
    urllib.request.urlretrieve(url, filename)

    for line in open(filename,encoding="utf-8").readlines():
        if "ネットワーク障害について" in line:
            datetimes.append(datetime.datetime.strptime(line.replace(" ","")[6:17], '%Y年%m月%d日'))
print(len(datetimes))

横軸に発生年月日、縦軸に障害累積回数を取って、どの時期に障害が多いのか確認する。

import matplotlib.pyplot as plt
plt.plot(list(reversed(datetimes)),[i+1 for i in range(len(datetimes))])
plt.xlabel("発生年月日")
plt.ylabel("障害累積回数")
plt.grid()
plt.show()

2018年は障害発生件数は比較的少ない他は、定常的に障害は発生しているようだ。

月ごとに発生回数のヒストグラムを作って、気温が高い夏に障害が多いのか、気温が低い冬に障害が少ないのか確認する。

plt.hist([dt.month for dt in datetimes], bins=12,range=[0.5,12.5])
plt.xticks([i for i in range(1,13)],[f"{i}月" for i in range(1,13)])
plt.xlabel("発生月")
plt.ylabel("発生回数")
plt.show()

2月は有意に少ない気がするが、他の月で傾向を見出すことはできなかった。

障害の内容が多種多様で、単なるハングアップから電源故障まで含むため、電源故障などの機器関連のみ抜き出せば、もう少し季節変動は見えたかもしれない。

障害発生から回復までの時間を取得する

import urllib.request
import datetime

suffix_urls = [""]
for i in range(2,14):
    suffix_urls.append(f"page{i}/")
datetimes = []
months = []
for suffix_url in suffix_urls:
    url = f"https://icts.nagoya-u.ac.jp/ja/information/trouble/{suffix_url}"
    filename = "a.html"
    urllib.request.urlretrieve(url, filename)

    for line in open(filename,encoding="utf-8").readlines():
        if "ネットワーク障害について" in line:
            url2 = "https://icts.nagoya-u.ac.jp"+line.split("\"")[1]
            filename2 = "b.html"
            urllib.request.urlretrieve(url2, filename2)
            lines2 = open(filename2,encoding="utf-8").readlines()
            for j, line2 in enumerate(lines2):
                if "障害発生時間" in line2:
                    if "~" not in lines2[j+1]:
                        print(lines2[j+2].replace("</p>",""),end="")
                    else:
                        print(lines2[j+1].replace("</p>",""),end="")
print(len(datetimes))

出力結果はここには載せない。気になる人は各自実行してみて欲しい。

Python + NumPy: 乱数生成のシード(RandomState)をちゃんと管理する

qiita.com

NumPyでは、発生する乱数をシードを使ってコントロールすることが可能で、グローバルシードを使う方法と、乱数状態を保持したクラスを使う方法がある。

グローバルシードを変更する np.random.seed(seed) と、乱数状態クラスを取得するための np.random.RandomState(seed) で、検証用のコードを書いてみた。なお、乱数の現在の状態を得るときは get_state() 関数を使う。

import numpy as np

seed=2

np.random.seed(seed)
print("グローバルシードを初期化 seed =",seed)

state = np.random.get_state()
print("状態取得",state[0],state[1][0],state[2:])

print("乱数発生",np.random.rand())

state = np.random.get_state()
print("状態取得",state[0],state[1][0],state[2:])

print()

rs = np.random.RandomState(seed)
print("ローカルシードを初期化 seed =",seed)

state = rs.get_state()
print("状態取得",state[0],state[1][0],state[2:])

print("乱数発生",rs.rand())

state = rs.get_state()
print("状態取得",state[0],state[1][0],state[2:])

どの環境で実行しても、

グローバルシードを初期化 seed = 2
状態取得 MT19937 2 (624, 0, 0.0)
乱数発生 0.43599490214200376
状態取得 MT19937 1091773355 (2, 0, 0.0)

ローカルシードを初期化 seed = 2
状態取得 MT19937 2 (624, 0, 0.0)
乱数発生 0.43599490214200376
状態取得 MT19937 1091773355 (2, 0, 0.0)

と出力されるはずである。

なお、seed = 2 で初期化した直後の state[1] を全て書くと

[         2 3624866507  846688490 2733819477 1447939927 2751582963
  523383323 3998955438 1552125001  468770929 3174827423  721042428
 3535242872 3938247060 2422198177 2908759646 3086638684 3470362919
 4130561862 1461945164 1969357429 3290893273 2478135064 2883471193
  442423807 3020847284 3712550376 2942769362 3450083884 2131393448
 1760508619 1079921617 2926983472 3248378587 2176075610 3115728347
  149850049 1325979722 2018146493 3978444371 1443587512 3757974310
 4061506243 3324845547  103150772 1574559537 4189437470  909936288
 4165778768  954014960 3460881122 1208468472 4206858865 2349975855
 1139307511  612263493 4068230769  836291891 3502129241 3858075069
  720214578 1721800695 1583913548 3992130976 1391510671 3658451527
 3506172182 3137165964 3175921738  668679341 2981266823 1736825536
 3777468269 2427100335 2783143307 3066847320  413674958 3816771731
 3611810078 2387881408 1275125722 4146169016 4116607257 1753188757
 2919070648  642669239 3128342921 3505362734 3299392025 4041168539
  194806354 2126204085 1155372896 3507921570 4183642851  390828735
  965876923 2502060328 3640719348  500564950 2958494162  726085493
 3551417999 1038946723 2355966135 2266502354 2169734778  103957955
 2239750491 3497757066 3960902779 3919891399 2708228292  672443535
 3966344157 3039186185 3064219851 3275966658 1958870427 3363207225
 1051359578 1034488571 4218333313 3731099077  826434970 3066527807
  503397775 2203570410 1724486664 1174354190  379816557  949285508
 1512592024 4098586594 4162106699 1244073455 4049526126  149401994
 2486298108 2645994433 4048691067 2760221413 1429726385  155150847
  263670059 1839522184 1971211423 2988886249 2906654795 2763085410
 3184555126 2428241499 1912808117 2344289693 1914982229  303251647
 3363291895 1933606625 4101500670 3877744496 2270413823 1339668594
 2478854401 2974781906 1176295348 1522792462 3413763473 1643164993
 1773585128  600348566 3831305432 3739184658 2927920737  921239740
 2960516826  653634535 4282151379 3899733953 2848969276  506667305
 2064254689 2562889493 1027730377 2964297732 1573662998 1400221132
  604719003 2600612834 1086412828 2839169326 2237360922 2533955895
 2173128617 3585208888  225649929 3162402384 2737400094  945452241
   32174139 4057815822 3681149417  538074139 1665518193 1290718971
 1204908910 1752298136 1576507435 1950581601  395392432 3117768513
 1610510657  582872851 3882188371 3486754405 3612818708   79522282
 3654222890  554014214 2795988792 2699019453  768333623 1041848464
 1012203950 3207710853 2191875363 2180736998 3657761238 2596779756
 1513261770  298083836  921966930 2568544321  480859735 4268007484
  757985477  442264484 4253891744  178333756 3653303450 2759887436
 1625051062 2140496676 2583817099  109920000 3849378804  587135080
 3965692414 3584128712 2496889359 4066788634 3439502551 1041005727
 2784445111 4231212902 3427909847  306088099  708846671 1357224748
   81546947 1013346034 3328657534 1140162390 1617022297  451026111
 4191194979 2544038377 3358406081 2078869653 1150523248 3663542690
 2315642515  380539204  125351140 2913559301 3540242901 1353142401
 2409055124 3136483651 3384539323 2622200239 1322006361   60288721
 1836454543 3378860321 1028797830  784192763 3200015141 4280866434
 3480136453 1163505791 1464340696   39962560 3132723172 4035042531
 3968861574 2859031456 1997068818 2914465192 3299223356  658977798
 2383986058 3689451477 1855211164  613736992 3399505872  721916528
 1214953314 1978881602 1644160163  868083487 1886701425 2430443623
 3238614545  440178771 4122051833  168378845 4075897453  784209059
 2122211469 1029375099 3987445959 1993353877 3774218662  259191900
  463709584 3256509205 2970214132   32546389 3103585745  404744328
 1405060850 3194833450 2229644052 1287997179 1168559088 4253642852
  456227059  115590192  959060034 2523194717 2921813711 2660496694
 1661244122 3949077694  546088169 3927753286  758289299 4076004698
 2818580089  297327332  391012690 1698646201 1827285496 2057818526
   45927453 4001566164 1358792247  340974259 3667858693 2116911813
 1716223164 3870579962  565768359 1330729550 1682620823  640788379
 2642810005  722202626 3342124090 1046670574 1884608856  926507664
 4048745540 1013276792 1474775246 3319759394  683213181 3644095690
 4261137863 1171395535 1797066690 2184942700 3136450276 2403915837
   19925339 3530712424  693454777 4095102592 2791708211 2407398362
  191034814 1405068157 1635573108 3002727858 2503231482 1877745251
 2971168822 4095482385 2810703016  116218529  663481365 3424543450
  181543471 2088477982 3758013903 4087811857 2306315696 1797804753
 3739876776 1052678416 2961962218  530565411 1084430187  600176751
 2920770665 1395657430 2275563379 1945041462 2897291349 3202120438
 2541033448 3395594999 1440745450  554493022 3604900030 1693940026
 2048944369 2776750171 1064369609 1004778234 3631442000  111252846
 2495206422 1549492629 1917957398 2076623558 2065193271 3503495683
  305349558 1064312453 4180949297 3969143923 1187332586 1837549170
 1199888411  801451519  878812761 3366984668 1132353467  467740707
 1531220625 1986847891 3611709790   97771638 2599466580  921040821
 3695024433  598160515 1067897977 1070962312 4203425140 3330187456
 3176319421  400150826 1778712930 1332512992  356084631 3210494822
 3878886472 4217065068   17324705 1194519644  169871241  594168806
 4141353368 2812067586 1007574236 3143655593 1231988821  990469891
  731912463 3875269580 3107901837 3268884046   51553093 3265091870
 2827988823  215460720 1810208024  235746758 2640785160 3542608605
 3212999042 3019862381 2776967609 2566622134 3136123892 3329714175
 1537598558 1111683694 1831866047 1527480555 1374492744 2861082564
 2283612438 3725126365 1934064016  976785712 1985214956 1901883774
 3769548825 1651881537 3535350080 2666996848  287335676  716617071
 1798766031 2468831819 2653018579 2711800444 2500416446 1688315701
 3098398094 2389050951  259942469 3096384070 2434046690    2292847
 2179858651  860453038 1112784056  806844432 2431840356 3444142163
  491004326 4117569173 1468356934 1230256668 3558831243  360563651
 1723242763 2771827983  555715903 1690465018 3642536487  312708693
 2718279339 3439997904 3381092195  269681413 3040927775 2538795672
 2982421226 3983011761 4012234084 3807531726 4152713741 2675650483
 3522595843 2603695663 1139711729  415442657  821281271 1153021862
 4196324887 2034707481 1146437550  645908802 3958160194 1331866590
 3875612725 2742273417 2083607827 1367028823 3405229612  926409162
  359084274  428200891 1580768777  723559275  390184571 1966930892
 1114312487 4049752645 2713447910  313516349 1670939739 2155242445
  415771127 1060028352  810510606 1441723093 3674555380 1244922308
 1138369035 2811021893 4005351255 1455276665 3073897390 1949184307
 4276825362 3579699982  501713787 1784405218 3987251947 1460563941
 1372531538 1479574814 2050694811 2867429155 3373192039  968223959
 1315340855  563759027 1277302533 2564831995 4174743461  794844391
 1091389581 3521416615  779397152 3672078093  440563956  765480627]

となる。メルセンヌ・ツイスターのページ数は 19937bit なので、32bit整数にして624個のメモリ空間を使うことになる。