# sudo apt install python3-pip
# python3 -m venv venv
# cd venv
# . bin/activate
# pip3 install numpy matplotlib
import numpy as np
import matplotlib.pyplot as plt
n = 100
t = np.linspace(0, 1.0, n)
y = np.zeros(n)
y[0] = 5 # initial boundry condition
# x(t) := a*exp(w*t) - b => b = a*exp(0) - x(0), dy/dt = w*a*exp(w*t) = w * [x(t) + b]
w = -5 # weight, - to decay, + to explode
a = 10 # amplitude
b = a - y[0] # bias
x = a * np.exp(w*t) - b # time evolution function, x(0) = a - b => b = a - x[0], 常微分 dx/dt = w*(x +b)
df_dt = lambda t, x: w * (x + b) # def df_dt(t, x): return w * (x + b)
z = y
dt = t[1] - t[0]
for i in range(n - 1): # to update [i + 1]
# 1. 4th Order Runge-Kutta method to solve y for dy/dt = f(y, t)
k1 = df_dt(t[i], y[i])
k2 = df_dt(t[i] + dt / 2, y[i] + dt * k1 / 2)
k3 = df_dt(t[i] + dt / 2, y[i] + dt * k2 / 2)
k4 = df_dt(t[i] + dt , y[i] + dt * k3 )
y[i + 1] = y[i] + dt * (k1 + 2 * k2 + 2 * k3 + k4) / 6
# 2. 1st Order Euler method to solve z for dz/dt = f(z, t)
z[i + 1] = z[i] + dt * df_dt(t[i], z[i])
fig = plt.figure()
plt.subplot(311)
plt.plot(y)
plt.ylabel('y: RK4')
plt.subplot(312)
plt.plot(z)
plt.ylabel('z', rotation=75)
plt.yticks(rotation=90)
plt.subplot(313)
plt.plot(x)
plt.ylabel('x:=a*exp(wt)-b')
plt.xlabel(f"t, a={a}, b={b}, w={w}")
plt.show()
2026年7月31日 星期五
用 python 解簡單的常微分方程式
2026年2月24日 星期二
用 emcc compiler 寫簡單的 openGL 繪圖, 讓瀏覽器也能觀看
在瀏覽開啟 html 檔, 裡面除了可以用 javascript 語言來運行openGL ES, 也能用 wasm 語言來運作 , 透過 emcc 編譯器可以將 c 語言翻譯成 wasm, 安裝方式詳如 Emscripten 官網:
https://emscripten.org/docs/getting_started/downloads.html
在 linux 系統上, 我將 emsdk 安裝到 ~/project/emsdk 目錄內, 底下是簡單的 Makefile 用來將檔案輸出到 ramdisk (/dev/shm)內, 只要執行 make 就能產生包含用 g++ 編譯的可執行檔 /dev/shm/main.out, 加上可以讓瀏覽器開啟的網頁(/dev/shm/main.htm) 和 javascript 執行檔 (/dev/shm/main_wasm.js)
#Makefile
# := 變數指定一次
# = 變數可以重複指定
# 所有來源 : $^
# 第1來源 : $<
# 目標主名+副名: $@
# 目標主名 : $*
c_SRC := gltest.cpp
js_HTML := main.htm
js_WASM := main_wasm.js
path_DST:= /dev/shm
path_SRC:= $(shell pwd)
cc_LIBs := -I include -l GL -l glut
em_LIBs := -I include -s WASM=1 -s LEGACY_GL_EMULATION=1 -s USE_WEBGL2=1 -s SINGLE_FILE -s USE_FREETYPE=1
rd_HTML := $(path_DST)/$(js_HTML)
rd_WASM := $(path_DST)/$(js_WASM)
rd_EXE := $(path_DST)/main.out
define html_content
<!DOCTYPE html>
<html><head><meta charset="utf-8"><meta http-equiv="Content-Type" content="text/html; charset=utf-8"></head><body>
<center>
<canvas id="canvas" oncontextmenu="event.preventDefault()"></canvas>
<script type="text/javascript">var id_from_canvas = document.getElementById("canvas");var Module = {canvas: id_from_canvas};</script>
<script src="$(js_WASM)"></script>
</center>
</body></html>
endef
all: $(rd_HTML) $(rd_WASM) $(rd_EXE)
@echo open $< to run $(rd_WASM) in browser
run: $(rd_EXE)
$(rd_EXE)
$(rd_WASM): $(c_SRC)
cd ~/project/emsdk && . emsdk_env.sh && cd $(path_SRC) && emcc $^ $(em_LIBs) -o $@ && echo " "
$(rd_HTML):
$(file > $@, $(html_content))
$(rd_EXE):$(c_SRC)
@g++ $^ $(cc_LIBs) -lfreetype -o $@
clean:
rm -f $(rd_WASM) $(rd_HTML) $(rd_EXE)
簡單用 c++ 寫一個繪圖程式 :
// gltest.cpp:
#include <stdio.h>
#include <string.h>
#include <stdlib.h>
#include <math.h>
#include <unistd.h>
#include <GL/glut.h>
struct ColorRGB { float r,g,b; };
const ColorRGB colors[] = {
{.r=1, .g=1, .b=1},
{.r=1, .g=1, .b=0},
{.r=1, .g=0, .b=1},
{.r=1, .g=0, .b=0},
{.r=0, .g=1, .b=1},
{.r=0, .g=1, .b=0},
{.r=0, .g=0, .b=1},
{.r=0.5, .g=0.5, .b=0.5}
};
int size_n = sizeof(colors)/sizeof(colors[0]);
void draw_circle (float cx, float cy, float radius, ColorRGB c = colors[0], int max_segments = 32) {
const double d_theta = M_PI * 2 / max_segments;
double theta = 0;// initial θ
int segments = max_segments;// lines to draw
glColor3f(c.r, c.g, c.b);
glBegin(GL_LINE_LOOP);// GL_TRIANGLE_FAN or GL_LINE_LOOP, to close loop
while (segments -- > 0) {
glVertex2f(cx + radius * cos(theta), cy + radius * sin(theta));
theta += d_theta;
}
glEnd();
}
void draw_line(float x0, float y0, float x1, float y1, ColorRGB c = colors[0]){
glColor3f(c.r, c.g, c.b);
glBegin(GL_LINES) ; // to draw one line
glVertex2f(x0, y0); // first point
glVertex2f(x1, y1); // second point
glEnd(); // end drawing
}
void update_loop(int parameter){
printf("parameter = %8d\n", parameter);
glutPostRedisplay();// send event to redraw
glutTimerFunc(1000, update_loop, parameter + 1);// continue to run update_loop again after 1 second
}
int main(int argc, char** argv) {
glutInit(&argc, argv); // Initialize GLUT
glutInitDisplayMode(GLUT_SINGLE | GLUT_RGB); // Set display mode
glutInitWindowSize(800, 600); // Set window size
glutCreateWindow("GLUT Example"); // Create window
atexit([]() {
printf("size_n = %-8d, bye~bye.\n", size_n);
});
glutDisplayFunc([]() {// when redraw event happens
static int current = 0;
current = (current + rand()) % size_n;
glClear(GL_COLOR_BUFFER_BIT);
glColor3f(0, 0, 0);
draw_circle( 0, 0, 0.5, colors[current]);
draw_line (-1, 1, 1, -1, colors[(current + 1) % size_n]);
draw_line (-1, -1, 1, 1, colors[(current + 2) % size_n]);
glFlush(); // Flush drawing command buffer. If using double buffering (GLUT_DOUBLE), use glutSwapBuffers();
}); // Register display callback
update_loop(0); // begin to send redraw event
glutMainLoop(); // Enter GLUT event processing loop
return 0;
}
後記: 如果 make 運行時出現錯誤, 有可能是 Makefile 內執行命令前面的縮排字元(\t: Tab 按鍵)被空白字元(' ': Space 按鍵)取代了, 只要用編輯器將它修正回縮排字元就能正常運作了.
2025年12月24日 星期三
在 linux 系統下簡單的 tar 檔案讀寫程式
參考網站: https://github.com/calccrypto/tar/tree/master,
改寫成我想用的: listtar.cpp
#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include <time.h>
#include <fcntl.h>
#include <errno.h>
#include <unistd.h>
#include <sys/stat.h>
#include <sys/types.h>
#include <grp.h>
#include <pwd.h>
#include <dirent.h>
#define debug_printf(fmt, ...) fprintf(stderr, fmt, ##__VA_ARGS__)
typedef struct Link_list_meta TarLinkList;
struct Link_list_meta {
union {
char block[512];// metadata
union {
struct {// Pre-POSIX.1-1988 format
char name[100]; // file name
char mode[8]; // permissions
char uid[8]; // user id (octal)
char gid[8]; // group id (octal)
char size[12]; // size (octal)
char mtime[12]; // modification time (octal)
char check[8]; // checksum of the block, with spaces in the check field while calculation is done (octal)
char link; // link indicator
char link_name[100]; // name of linked file
};
struct {// UStar: Unix Standard TAR format (POSIX IEEE P1003.1)
char old[156]; // first 156 octets of Pre-POSIX.1-1988 format
char filetype; // file type
char also_link_name[100]; // name of linked file
char ustar[6]; // ustar\0
char version[2]; // #Version
char owner[32]; // user name (string)
char group[32]; // group name (string)
char major[8]; // device major number
char minor[8]; // device minor number
char prefix[155];
};
};
};
TarLinkList *next;
ssize_t append(int fd, void *buf, int n) { return write(fd, buf, n); }// todo: append buf into fd at the end
void block_update(int fd, char *filename, off_t filesize, mode_t filemode = 0, time_t *ct = nullptr) {
if (fd < 0) return;
memset(check, ' ', sizeof(check));// init string, It must be empty before caculation.
sprintf(name, "%s" , filename);
sprintf(mode, "%07o", filemode > 0 ? filemode & 0777 : 0664);
if (filesize > 0) sprintf(size , "%011o",(unsigned int)filesize);
else memset(size , '0', sizeof(size));
if (ct) sprintf(mtime, "%011o",(unsigned int)*ct);
else { // using current time if not provide.
time_t now;
time(&now);
sprintf(mtime, "%011o",(unsigned int)now);
}
int n = sizeof(block), checksum = 0;
for (int i = 0; i < n; i++) checksum += (unsigned char)block[i];// caculate checksum in the block
sprintf(check, "%07o", checksum);
append(fd, block, n);
}
int open_ram2tar(char *dir_name, char *create_name=nullptr) {//todo: validate dir_name
char backup_name[strlen(dir_name) + 16];
if (create_name == nullptr) {
sprintf(backup_name, "_%s.tar", dir_name);
create_name = backup_name;
}
int tar_fd = open(create_name, O_RDWR | O_TRUNC | O_CREAT, S_IRUSR | S_IWUSR);
if (tar_fd < 0) return -1;
block_update(tar_fd, dir_name, 0, 0775);// chmod ug+rwx o+rx dir_name
filetype = '0';// change to normal file.
debug_printf("createe file: %s\n", create_name);
return tar_fd;
}
Link_list_meta(bool is_directory = true) {// constructor to initialize the block
memset(block, 0, sizeof(block));
version[0] = '0', version[1] = '0';
filetype = is_directory ? '5' : '0';// todo: other type
uid_t user_id = getuid();
sprintf(uid, "%07o", user_id);
sprintf(gid, "%07o", getgid());
sprintf(ustar, "%s", "ustar");
struct passwd *pwd = getpwuid(user_id);// user info get from UID
if (pwd) {
struct group *grp = getgrgid(pwd->pw_gid);
sprintf(owner, "%s", pwd->pw_name);
sprintf(group, "%s", grp ? grp->gr_name: "None");
}
}
};
bool is_empty(char *buffer, int n) { // make sure first n's data in buffer are all 0s
for (int i = 0; i < n; i ++) if (*buffer ++) return false;
return true;
}
long int o2l(char *octal_str, int n) {// 8 進位轉長整數, todo: negative number
long int val_long = 0l;
for (int i = 0; i < n; i ++, octal_str ++) {
if (*octal_str == 0) break;
val_long <<= 3;
val_long |= *octal_str - '0';
}
return val_long;
}
void dir2tar(const char *foldername, char *create_name = nullptr) {
char *path2folder = (char *)foldername;
if (*path2folder == '.') {
path2folder ++;
if (*path2folder == '.') path2folder ++;
}
if (*path2folder == '/' ) path2folder ++;
else path2folder = (char *)foldername;
while (*path2folder == '/') path2folder ++;// remove another '/'
struct stat file_st;
if (lstat(path2folder, &file_st)!=0 || (file_st.st_mode & S_IFMT)!=S_IFDIR) return;
DIR *cd = opendir(path2folder);// change into the directory
if (!cd) return; // make sure user has rights to access
TarLinkList metadata;
int tar_fd = metadata.open_ram2tar(path2folder);
if (tar_fd > 0) {
char fd_buf[512];
struct dirent *temp;
while ((temp = readdir(cd))) { // todo: to proceed child directory
if (temp->d_name[0] == '.') continue;// skip . and ..
char temp_fullname[strlen(path2folder) + strlen(temp->d_name) + 2];// + '/' and EOS
sprintf(temp_fullname, "%s/%s", path2folder, temp->d_name);// fullname
lstat(temp_fullname, &file_st);
if ((file_st.st_mode & S_IFMT) != S_IFREG || file_st.st_size <= 0) continue; // todo: other type support
int temp_fd = open(temp_fullname, O_RDONLY);
if (temp_fd < 0) continue;
off_t temp_len = file_st.st_size;// to append file content
off_t zeros_pad = temp_len % 512;// check remain
if (zeros_pad) zeros_pad = 512 - zeros_pad;// number of zeros need to pad
metadata.block_update(tar_fd, temp_fullname, temp_len, file_st.st_mode, &file_st.st_mtim.tv_sec);
debug_printf("%s: size = %ld, checksum = %6s\n", temp_fullname, temp_len, metadata.check);
while (temp_len > 0) {
ssize_t l = read(temp_fd, fd_buf, temp_len > 512 ? 512 : temp_len);
if (l <= 0) break;//todo: error correction
metadata.append(tar_fd, fd_buf, l);
temp_len -= l;
}
close(temp_fd);
if (zeros_pad) {
memset(fd_buf, 0, zeros_pad);
metadata.append(tar_fd, fd_buf, zeros_pad);
}
}
memset(fd_buf, 0, 512);// need 2 block of zeros in the end for tar file
for (int i = 0; i < 2; i ++) metadata.append(tar_fd, fd_buf, 512);
close(tar_fd);
}
closedir(cd);
}
void list_tarfile(const char *tar) {
int fd = open(tar, O_RDONLY);
if (fd > 0) {// in linux: stdin = 0, stdout = 1, stderr = 2
TarLinkList *archive = nullptr;// start
TarLinkList **tarlist = &archive;// get pointer of archive
int block_size = sizeof(archive->block);
while (true) {
TarLinkList *temp = (TarLinkList *)calloc(1, sizeof(TarLinkList));// 分配空間並初始為 0
if (temp == nullptr) break;
if (read(fd, temp->block, block_size) != block_size) {// to read 512 bytes metadata
debug_printf("讀取錯誤,忽略!\n");
free(temp);
break;
}
if (is_empty(temp->block, block_size)) {// EOF, enough to stop
if (read(fd, temp->block, block_size) == block_size) {// check 2nd EOF
if (is_empty(temp->block, block_size)) {
debug_printf("正常檔尾,結束:\n");
}
}
free(temp);
break;
}
*tarlist = temp;// fill temp as current entry
tarlist = &temp->next;// to fill for next time
long int goahead = o2l(temp->size, 11);// 檔案長度: 8 進位(12 bytes)
int r = goahead % 512; // 取餘數
if (r) goahead += 512 - r;// 若非 512 倍數, 無條件補滿成 512 倍數
if (lseek(fd, goahead, SEEK_CUR) < 0) { // 前進到下個位置
debug_printf("前進錯誤,忽略!\n");
break;
}
}
*tarlist = nullptr;// end of List
while (archive) {// list and free
time_t t = o2l(archive->mtime, 11);// 更新時間
struct tm *ct = localtime(&t);
debug_printf("%s@%s\t", archive->owner, archive->group); // 使用者@群組
debug_printf("%d-%02d-%02d:%02d.%02d\t",
ct->tm_year + 1900,
ct->tm_mon + 1,
ct->tm_mday,
ct->tm_hour,
ct->tm_min
);// 年-月-日-時:分
switch (archive->filetype) {
case '0':
debug_printf("%ld (bytes)", o2l(archive->size, 11));// 檔案長度
break;
case '1': case '2':
debug_printf("檔案連結");
break;
case '3': case '4':
debug_printf("裝置檔案-%04ld::%04ld-", o2l(archive->major, 7), o2l(archive->minor, 7));// 設備編號
break;
case '5':
debug_printf(" <目錄> ");
break;
case '6':
debug_printf("先進先出");
break;
default:
debug_printf("????");
break;
}
debug_printf("\t<- (%6s) %-32s\n", archive->check, archive->name);// 檔名
TarLinkList *temp = archive -> next;// remove later
free(archive);
archive = temp;
}
close(fd);
}
}
void dump_tarfile(const char *tar, const char *filename) {
int fd = open(tar, O_RDONLY);
if (fd > 0) {
char fd_buf[512];// as buffer
TarLinkList *archive = (TarLinkList *)fd_buf; // point to fd_buf
long int goahead = 0l;
while (lseek(fd, goahead, SEEK_CUR)>= 0 && read(fd, fd_buf, 512) == 512 && !is_empty(fd_buf, 512)) {
if (strcmp(archive->name, filename) == 0) {
long int filesize = o2l(archive->size, 11);
while (filesize > 0) {
int l = read(fd, fd_buf, (filesize > 512) ? 512 : filesize);
if (l <= 0) continue;
for (int i = 0; i < l; i++) debug_printf("%c", fd_buf[i]);
filesize -= l;
}
break;
}
goahead = o2l(archive->size, 11);// 檔案長度: 8 進位(12 bytes)
int r = goahead % 512; // 取餘數
if (r) goahead += 512 - r;// 若非 512 倍數, 無條件補滿成 512 倍數
}
close(fd);
}
}
int main(int argc, char *argv[]) {
if (argc > 1 && argv[1]) {
struct stat file_st;
if (lstat(argv[1], &file_st) == 0) {// make sure argv[1] file exists.
if ((file_st.st_mode & S_IFMT) == S_IFDIR) {
dir2tar(argv[1]);// createe tar file to store all files in argv[1] which is a directory.
} else {// todo: make sure argv[1] is a tar file
if (argc > 2 && argv[2]) {
printf("===%s:%s===\n", argv[1], argv[2]);
dump_tarfile(argv[1], argv[2]); // to dump argv[2] in argv[1]
printf("\n=== EOF ===\n");// end of file
} else {
list_tarfile(argv[1]);// list all files in tar
}
}
}
}
return 0;
}
一個將 ram 資料寫入 tar file 測試程式: test_ram2tar.c
#include <stdio.h>
#include <string.h>
#include <time.h>
#include <fcntl.h>
#include <unistd.h>
#include <grp.h>
#include <pwd.h>
#define debug_printf(fmt, ...) fprintf(stderr, fmt, ##__VA_ARGS__)
typedef struct Link_list_meta TarLinkList;
struct Link_list_meta {
union {
char block[512];// metadata
union {
struct {// Pre-POSIX.1-1988 format
char name[100]; // file name
char mode[8]; // permissions
char uid[8]; // user id (octal)
char gid[8]; // group id (octal)
char size[12]; // size (octal)
char mtime[12]; // modification time (octal)
char check[8]; // checksum of the block, with spaces in the check field while calculation is done (octal)
char link; // link indicator
char link_name[100]; // name of linked file
};
struct {// UStar: Unix Standard TAR format (POSIX IEEE P1003.1)
char old[156]; // first 156 octets of Pre-POSIX.1-1988 format
char filetype; // file type
char also_link_name[100]; // name of linked file
char ustar[6]; // ustar\0
char version[2]; // #Version
char owner[32]; // user name (string)
char group[32]; // group name (string)
char major[8]; // device major number
char minor[8]; // device minor number
char prefix[155];
};
};
};
TarLinkList *next;
ssize_t append(int fd, void *buf, int n) { return write(fd, buf, n); }// todo: append buf into fd at the end
void block_update(int fd, char *filename, off_t filesize, mode_t filemode = 0, time_t *ct = nullptr) {
if (fd < 0) return;
memset(check, ' ', sizeof(check));// init string, It must be empty before caculation.
sprintf(name, "%s" , filename);
sprintf(mode, "%07o", filemode > 0 ? filemode & 0777 : 0664);
if (filesize > 0) sprintf(size , "%011o",(unsigned int)filesize);
else memset(size , '0', sizeof(size));
if (ct) sprintf(mtime, "%011o",(unsigned int)*ct);
else { // using current time if not provide.
time_t now;
time(&now);
sprintf(mtime, "%011o",(unsigned int)now);
}
int n = sizeof(block), checksum = 0;
for (int i = 0; i < n; i++) checksum += (unsigned char)block[i];// caculate checksum in the block
sprintf(check, "%07o", checksum);
append(fd, block, n);
}
int open_ram2tar(char *dir_name, char *create_name=nullptr) {//todo: validate dir_name
char backup_name[strlen(dir_name) + 16];
if (create_name == nullptr) {
sprintf(backup_name, "_%s.tar", dir_name);
create_name = backup_name;
}
int tar_fd = open(create_name, O_RDWR | O_TRUNC | O_CREAT, S_IRUSR | S_IWUSR);
if (tar_fd < 0) return -1;
block_update(tar_fd, dir_name, 0, 0775);// chmod ug+rwx o+rx dir_name
filetype = '0';// change to normal file.
debug_printf("createe file: %s\n", create_name);
return tar_fd;
}
Link_list_meta(bool is_directory = true) {// constructor to initialize the block
memset(block, 0, sizeof(block));
version[0] = '0', version[1] = '0';
filetype = is_directory ? '5' : '0';// todo: other type
uid_t user_id = getuid();
sprintf(uid, "%07o", user_id);
sprintf(gid, "%07o", getgid());
sprintf(ustar, "%s", "ustar");
struct passwd *pwd = getpwuid(user_id);// user info get from UID
if (pwd) {
struct group *grp = getgrgid(pwd->pw_gid);
sprintf(owner, "%s", pwd->pw_name);
sprintf(group, "%s", grp ? grp->gr_name: "None");
}
}
};
int main(int argc, char *argv[]) {
TarLinkList metadata;
char dir_name[16] = {"bb"};
int tar_fd = metadata.open_ram2tar(dir_name);// create folder first
if (tar_fd > 0) {
char fd_buf[512], ram_name[100];
for(int i = 0; i < 10; i ++) { // to create 10 example files
sprintf(ram_name, "%s/%d", dir_name, i);// combine folder name with specific name as ram file name
sprintf(fd_buf, "%s:%d", ram_name, i + 1);// fill content for the ram file
metadata.block_update(tar_fd, ram_name, strlen(fd_buf));// append metadata into tar
metadata.append(tar_fd, fd_buf, 512);// append 512 bytes of ram into tar
}
memset(fd_buf, 0, 512);
for (int i = 0; i < 2; i ++) metadata.append(tar_fd, fd_buf, 512);// append 2 block of zeros in the end
close(tar_fd);
}
return 0;
}
2025年11月29日 星期六
使用 linux 玩早期的 Turbo C
Turbo C 是早期在 dos (磁碟作業系統)底下的用來編譯 C 語言的編譯程式,可以上官網下載:
https://turbo-c.net/turbo-c-download/
為了讓它能在 linux 底下順利運作, 可以安裝 dosbox:
sudo apt install dosbox
或是上 dosbox 官網 https://sourceforge.net/projects/dosbox/files/dosbox/0.74-3/
下載原始程式, 自行編譯, 但要事先安裝必要的程式庫:
sudo apt install libsdl1.2-dev
解壓縮後, 只要在原始目錄底下運行
./configure && make
就會在 src/ 目錄下編譯出可執行檔(src/dosbox), 將它複製到任何需要在 dos下運作的程式目錄下, 伴隨 dosbox 可執行檔, 在 dosbox 所在目錄下, 可以自行編輯一個 dosbox.conf 將開機後要執行的命令放在裡面, 讓它自動執行開機後的執行命令, 省下許多打字的時間, 例如:
[autoexec]
mount c ~/project/dos/TURBOC3
path=c:\BIN
c:
dosbox 目前已經沒在更新, 若要使用仍在維護的 dosbox 版本, 另外有 dosbox-x, 可以上官網下載原始檔: https://github.com/joncampbell123/dosbox-x/releases
但要事先安裝許多必要的工具程式及程式庫:
sudo apt install automake nasm libncurses-dev libsdl-net1.2-dev libsdl2-net-dev libpcap-dev libslirp-dev fluidsynth libfluidsynth-dev libavformat-dev libavcodec-dev libavcodec-extra libswscale-dev libfreetype-dev libxkbfile-dev libxrandr-dev
解壓縮後, 只要在原始目錄底下運行
./build-debug
最後在 src/ 目錄下產生可執行檔(src/dosbox-x), 同 dosbox 可以自行編輯一個 dosbox.conf 將開機後要執行的命令放在裡面. 以後只要執行該目錄底下的 dosbox-x 就可
2025年11月7日 星期五
Samsung A16 移除雲端 app
adb pair #ip_address:#wifi_pair_tcp_port
adb connect #ip_address:#wifi_debug_tcp_port
adb shell pm list packages | grep "scloud"
adb shell pm uninstall --user 0 com.samsung.android.scloud
終於移除腦人的通知訊息, 真的很白目, 無言 ...
註:
1. #ip_address 是手機的 ip 位址
2. #wifi_pair_tcp_port 是要配對的 tcp 編號
3. #wifi_debug_tcp_port 是透過 WiFi 的 debugging tcp 編號
2025年9月11日 星期四
在 linux 上使用 qemu 玩 android
1. 安裝 qemu 及相關工具程式 : sudo apt install qemu-system-x86 qemu-utils
2. 上 andoid-x86 網站下載 android-x86_64-9.0-r2-k49.iso 檔 : https://sourceforge.net/projects/android-x86/files/
3. 事先建立好 20G 的虛擬機影像檔: qemu-img create -f qcow2 x86.qcow2 20G
4. 開機啟動 iso 檔, 需按照螢幕指示, 先創造並切割硬碟分割區, 最後將 android 系統安裝到虛擬機
qemu-system-x86_64 -enable-kvm -drive file=x86.qcow2,if=virtio \
-machine type=q35,vmport=off \
-display sdl,gl=on \
-audiodev pa,id=snd0 \
-device AC97,audiodev=snd0 \
-device virtio-vga-gl \
-device virtio-tablet \
-device virtio-keyboard \
-device qemu-xhci,id=xhci \
-net nic,model=virtio-net-pci \
-net user,hostfwd=tcp::4444-:5555 \
-cpu host -m 4096 -usb -smp 4 \
-cdrom android-x86_64-9.0-r2-k49.iso
5. 以後只要啟動虛擬機就可以了, 不再需要 iso 檔. 記得將 smp 數量降低, 避免全數 smp 被使用.
qemu-system-x86_64 -enable-kvm -drive file=x86.qcow2,if=virtio \
-machine type=q35,vmport=off \
-display sdl,gl=on \
-audiodev pa,id=snd0 \
-device AC97,audiodev=snd0 \
-device virtio-vga-gl \
-device virtio-tablet \
-device virtio-keyboard \
-device qemu-xhci,id=xhci \
-net nic,model=virtio-net-pci \
-net user,hostfwd=tcp::4444-:5555 \
-cpu host -m 4096 -usb -smp 2
2025年7月27日 星期日
使用 vscode 時, 改善滑鼠反應遲鈍的問題
按下 Manage 按鈕 Settings, 輸入 server, 儘量避免選項被啟用, 讓選項儘量 disable 或 off 例如:
Http: Proxy Strict SSL 不要勾選
C_Cpp: Code Folding 選擇 disable
C_Cpp: Suggest Snippets 不要勾選
2025年7月25日 星期五
簡單利用 sdl 載入 jpeg 檔, 描繪中文字, 線/圓繪圖
// sudo apt install libsdl2-dev libsdl2-image-dev libsdl2-ttf-dev
// g++ sdldraw.cpp -lSDL2 -lSDL2_image -lSDL2_ttf && ./a.out
// sdldraw.cpp
#include <unistd.h>
#include <SDL2/SDL.h>
#include <SDL2/SDL_image.h>
#include <SDL2/SDL_ttf.h>
void drawcircle(SDL_Renderer *renderer, float cx, float cy, float r) {
const int max_segments = 32;// large enough to smooth circle
const double d_theta = M_PI * 2 / max_segments;
double theta = d_theta;// 2nd θ
int px = cx + r;// 1st θ = 0
int py = cy;
int lines = max_segments - 1;// lines to draw
while (lines -- > 0) {
int nx = cx + r*cosf(theta);
int ny = cy + r*sinf(theta);
SDL_RenderDrawLine(renderer, px, py, nx , ny);
theta += d_theta;// next θ
px = nx;
py = ny;
}
SDL_RenderDrawLine(renderer, px, py, cx + r , cy);// close loop
}
int main(int argc, char** argv) {
if (SDL_Init(SDL_INIT_EVERYTHING) == 0) {
SDL_Window *xwin = SDL_CreateWindow("繪圖程式", 0, 0, 800, 600, SDL_WINDOW_RESIZABLE);
if (xwin) {
SDL_Renderer *renderer = SDL_CreateRenderer(xwin, -1, SDL_RENDERER_ACCELERATED);
if (renderer) {
SDL_Texture *bgPicture = IMG_LoadTexture(renderer, "snap.jpg");
SDL_RenderCopy(renderer, bgPicture, NULL, NULL);
SDL_RenderPresent(renderer);// show current image
SDL_SetRenderDrawColor(renderer, 255, 0, 0, 255);
const char *message = "準心";
const SDL_Color colorGreen = {.r=0, .g=255, .b=0, .a=255};
TTF_Font *ukai = (TTF_Init() == 0) ? TTF_OpenFont("./fonts/ukai.ttc", 32) : nullptr;
SDL_Rect target_cross;
SDL_Surface *msgSurface = ukai ? TTF_RenderUTF8_Blended(ukai, message, colorGreen) : nullptr;
SDL_Texture *msgTexture = SDL_CreateTextureFromSurface(renderer, msgSurface);
int &radius = target_cross.w; // alias
if (msgSurface) {
target_cross.w = msgSurface->w;
target_cross.h = msgSurface->h;
SDL_FreeSurface(msgSurface);
}
SDL_Event event;
while (true) { // event loop begin
usleep(1000);
SDL_PollEvent(&event);
if (event.type == SDL_QUIT) break;
if (event.type == SDL_MOUSEBUTTONDOWN) {
if (event.button.button == SDL_BUTTON_LEFT) {
SDL_RenderCopy(renderer, bgPicture, NULL, NULL);
int px = event.button.x;
int py = event.button.y;
if (msgTexture) {
target_cross.x = px - target_cross.w/2;
target_cross.y = py - target_cross.h/2;
drawcircle(renderer, px, py, radius);// green circle
SDL_RenderDrawLine(renderer, px, py - radius, px, py + radius);// red cross
SDL_RenderDrawLine(renderer, px - radius, py, px + radius, py);
SDL_RenderCopy(renderer, msgTexture, NULL, &target_cross);
}
SDL_RenderPresent(renderer);
printf("Left mouse is down @(%d,%d)\n", px, py);
}
} else if (event.type == SDL_WINDOWEVENT) {
if (event.window.event == SDL_WINDOWEVENT_RESIZED) {
SDL_RenderCopy(renderer, bgPicture, NULL, NULL);
SDL_RenderPresent(renderer);
}
}
}
if (bgPicture) SDL_DestroyTexture(bgPicture);
if (msgTexture) SDL_DestroyTexture(msgTexture);
if (ukai) TTF_CloseFont(ukai);
SDL_DestroyRenderer(renderer);
}
SDL_DestroyWindow(xwin);
TTF_Quit();
}
SDL_Quit();
}
return 0;
}
後記. 2025.07.29 改用 sdl3, 程式庫事先要從原始碼編譯並安裝, 上述程式修改並重新編譯:
// g++ sdl3draw.cpp -lSDL3 -lSDL3_image -lSDL3_ttf && ./a.out
// sdl3draw.cpp
#include <stdio.h>
#include <unistd.h>
#include <math.h>
#include <SDL3/SDL.h>
#include <SDL3_image/SDL_image.h>
#include <SDL3_ttf/SDL_ttf.h>
void drawcircle(SDL_Renderer *renderer, float cx, float cy, float r) {
const int max_segments = 32;
const double d_theta = M_PI * 2 / max_segments;
double theta = d_theta;
int px = cx + r;
int py = cy;
int lines = max_segments - 1;
while (lines -- > 0) {
int nx = cx + r*cosf(theta);
int ny = cy + r*sinf(theta);
SDL_RenderLine(renderer, px, py, nx , ny);
theta += d_theta;
px = nx;
py = ny;
}
SDL_RenderLine(renderer, px, py, cx + r , cy);// close loop
}
int main(int argc, char** argv) {
if (SDL_Init(SDL_INIT_EVENTS)) {
SDL_Window *xwin = SDL_CreateWindow("繪圖程式", 800, 600, SDL_WINDOW_RESIZABLE);
if (xwin) {
SDL_Renderer *renderer = SDL_CreateRenderer(xwin, nullptr);
if (renderer) {
SDL_Texture *bgPicture = IMG_LoadTexture(renderer, "snap.jpg");
SDL_RenderTexture(renderer, bgPicture, NULL, NULL);
SDL_RenderPresent(renderer);
SDL_SetRenderDrawColor(renderer, 255, 0, 0, 255);
const char *message = "準心";
const SDL_Color colorGreen = {.r=0, .g=255, .b=0, .a=255};
TTF_Font *ukai = (TTF_Init()) ? TTF_OpenFont("./fonts/ukai.ttc", 32) : nullptr;
SDL_FRect target_cross;
SDL_Surface *msgSurface = ukai ? TTF_RenderText_Blended(ukai, message, 0, colorGreen) : nullptr;
SDL_Texture *msgTexture = SDL_CreateTextureFromSurface(renderer, msgSurface);
float &radius = target_cross.w; // alias
if (msgSurface) {
target_cross.w = msgSurface->w;
target_cross.h = msgSurface->h;
SDL_DestroySurface(msgSurface);
}
SDL_Event event;
while (true) {
usleep(1000);
SDL_PollEvent(&event);
if (event.type == SDL_EVENT_QUIT) break;
if (event.type == SDL_EVENT_MOUSE_BUTTON_DOWN) {
if (event.button.button == SDL_BUTTON_LEFT) {
SDL_RenderTexture(renderer, bgPicture, NULL, NULL);
int px = event.button.x;
int py = event.button.y;
if (msgTexture) {
target_cross.x = px - target_cross.w/2;
target_cross.y = py - target_cross.h/2;
drawcircle(renderer, px, py, radius);// green circle
SDL_RenderLine(renderer, px, py - radius, px, py + radius);// red cross
SDL_RenderLine(renderer, px - radius, py, px + radius, py);
SDL_RenderTexture(renderer, msgTexture, NULL, &target_cross);
}
SDL_RenderPresent(renderer);
printf("Left mouse is down @(%d,%d)\n", px, py);
}
} else if (event.window.type == SDL_EVENT_WINDOW_RESIZED) {
SDL_RenderTexture(renderer, bgPicture, NULL, NULL);
SDL_RenderPresent(renderer);
}
}
if (bgPicture) SDL_DestroyTexture(bgPicture);
if (msgTexture) SDL_DestroyTexture(msgTexture);
if (ukai) TTF_CloseFont(ukai);
SDL_DestroyRenderer(renderer);
}
SDL_DestroyWindow(xwin);
TTF_Quit();
}
SDL_Quit();
}
return 0;
}
2025年4月3日 星期四
使用 python 簡單實現多 cpu 平行處理
import multiprocessing as mp
import time
def iso_task(k):
print(f"task {k} @{time.time()} sec: sleep for 1 second")
time.sleep(1)
print(f"task {k} @{time.time()} sec: finish.")
n = mp.cpu_count()
print(f"Total CPUs = {n}")
tasks = []
start_time = time.time()
for i in range(n): # prepare all tasks to run
task = mp.Process(target=iso_task, args=(i,))
tasks.append(task)
for i in range(n): # fire all tasks at the same time
tasks[i].start()
for i in range(n): # wait all tasks to finish
tasks[i].join()
dt = time.time() - start_time
print(f"{round(dt, 3)} sec elapsed")
2025年4月2日 星期三
使用 python 實現 chebyshev 多項式及內插法
不囉唆, 詳如以下代碼：
import numpy as np
def Cp(x, n): # 快速疊代法, 計算 1st kind Chebyshev 多項式
if (n == 0):
return 1
if (n == 1):
return x
pk_1 = 1
pn_1 = x
k = 1
while k < n :
pk = pn_1 * x * 2 - pk_1
pk_1 = pn_1
pn_1 = pk
k += 1
return pn_1
test_f = lambda x: np.exp(-x*x) # 測試函式
Ln = 20
px = [0.0] * Ln
py = [0.0] * Ln
theta = [0.0] * Ln
for k in range(Ln) :
theta[k] = np.pi * (k + 0.5) / Ln;# θₖ = np.pi * (k + 0.5) / Ln
px[k] = np.cos(theta[k]);# pxₖ = cos(θₖ)
py[k] = test_f(px[k]);# pyₖ = f(pxₖ) to be used in Chebyshev Interpolation
def Ci(i): # Coefficient, bind with pyₖ, θₖ, Ln
sum = 0
for k in range(Ln) :
sum += py[k] * np.cos(i * theta[k])
return sum * 2 / Ln # 2/n * Σₙf(pxₖ)*Tₖ(pxₖ), k = 0, 2, ... n - 1
def Chebyshev_interpolation(t): # Chebyshev Interpolation Polynomials
sum = Ci(0) / 2;# Σₙ Cₖ*Tₖ(x) - C₀/2 = Σₖ Cₖ*Tₖ(x) + C₀/2, k = 1, 2, ... n-1, C₀*T₀(x) = C₀
k = 1
while k < Ln :
sum += Ci(k) * Cp(t, k)
k += 1
return sum
for k in range(Ln) :
xk = px[k] + 0.1
c = test_f(xk) # 實際值
h = Chebyshev_interpolation(xk)# Chebyshev 合成值
print(f"x={xk}, 內插={h}, 實際值={c}, 誤差 = {h-c}")
2025年3月20日 星期四
用 python 實現定積分 ∫ₐᵇ f(t)dt
看了一些文章後, 自己手動寫了一些簡單的程式, 實現各種定積分的方式
import numpy as np
def Lp(x, n) : # evaluate order n-1 Legendre polynomials at x
if (n == 0) :
return 1.0
if (n == 1) :
return x
pn_1, pk_1 = x, 1.0 # 疊代初始化, 因為 n > 1, 所以至少疊代一次
k = 1 # 此刻從 k 開始, 直到 n - 1
while k < n :
pn_1, pk_1 = (pn_1*x*(2 * k + 1) - pk_1*k) / (k + 1), pn_1
k += 1 # # 因為 k=n-1 所以 n=k+1, 2*n-1 = 2*(k+1)-1 = 2*k+1
return pn_1 # (Lp(x, k)*x*(2*n-1) - Lp(x, k-1)*k) / n, k = n - 1
def d_Lp(x, n) : # derivative of order n-1 Legendre Polynomials at x
if (n == 0) :
return 0.0
if (n == 1) :
return 1.0;# LP′ₙ(x) = (−LPₙ(x)*x + LPₖ(x)) * n / (1−x*x), k=n-1
return (Lp(x, n - 1) - Lp(x, n) * x) * n / (1.0 - x*x)
# Multiple root finder algorithm for Legendre polynomial: https://publikacio.uni-eszterhazy.hu/3009/1/AMI_33_from3to13.pdf
def by_newton(f, a=1.0, b=2.0, n=10):
def newton_raphson(n, eps=1e-16): # using: xₖ = xₖ - f(xₖ) / [f'(xₖ) - f(xₖ) Σₖ 1/(xₖ - xᵢ)]
e = np.cos(np.array([np.pi*(k + 0.5)/n for k in range(n)])) #initial guess
for k in range(n) : # find all roots by newton raphson method
xk = e[k]
iteration = 0
while (iteration < 1000) :
iteration += 1
f = Lp(xk, n)
temp = f if (f > 0) else -f
if (temp < eps): # 收斂
break
sum_r = 0.0
for j in range(k) :# sum_r = Σₖ 1/(xₖ - xᵢ) to remove previouse root
delta = xk - e[j]
temp = delta if (delta < 0) else -delta
if (temp < eps) :
continue# skip singular value!
sum_r += 1.0 / delta
dx = f / (d_Lp(xk, n) - f * sum_r)
temp = dx if (dx > 0) else -dx
if (temp < eps) : # 收斂
break
xk -= dx # xₖ = xₖ - f(xₖ) / [f'(xₖ) - f(xₖ) Σₖ 1/(xₖ - xᵢ)]
e[k] = xk # final root update
xi = [e]
derivative = d_Lp(xi[0], n)
xi.append(2 / (derivative*derivative*(1 - xi[0]*xi[0])))
return np.array(xi)
xw = newton_raphson(n)
scale = (b - a) / 2.0 # linear transform tx: bias + scale * 1 = b => scale = (b - a) / 2
bias = (b + a) / 2.0 # linear transform tx: bias + scale *(-1) = a => bias = (b + a) / 2
tx = xw[0] * scale + bias # x -> tx
sum = f(tx).dot(xw[1])
return sum * scale
def by_jacobi(f, a=1.0, b=2.0, n=10):
def Jacobi_method(n, eps=1e-16): # using: xₖ = xₖ - f(xₖ) / [f'(xₖ) - f(xₖ) Σₖ 1/(xₖ - xᵢ)]
J = np.zeros((n, n))
d = len(J) - 1 # to fill into the jacobian tridiagonal matrix, trace = 0, and it is a symmetry matrix
for k in range(d): # fill
m = k + 1
J[k][m] = m / np.sqrt(4.0 * m * m - 1.0)
J[m][k] = J[k][m]
e = np.linalg.eigvals(J) # nxn jacobi matrix, solve eigenvalues
xi = [np.sort(e)] # eigenvalue is same as root of the order n-1 Legendre polynopmial
derivative = d_Lp(xi[0], n)
xi.append(2 / (derivative*derivative*(1 - xi[0]*xi[0])))# weight relative eigenvalue
return np.array(xi)
xw = Jacobi_method(n)
scale = (b - a) / 2.0 # linear transform tx: bias + scale * 1 = b => scale = (b - a) / 2
bias = (b + a) / 2.0 # linear transform tx: bias + scale *(-1) = a => bias = (b + a) / 2
tx = xw[0] * scale + bias # x -> tx
sum = f(tx).dot(xw[1])
return sum * scale
def by_trapezoid(f, a=1.0, b=2.0, n=50): # integral with the trapezoid method, 使用 梯形面積=(上底 + 下底)/2 積分
x = np.linspace(a, b, n) # total n pints include a, b
y = f(x) # total n
n -= 1 # split to n - 1 interval
delta_x = (b - a) / n
sum = y[1:n].sum() + (y[0] + y[n]) / 2.0
return sum * delta_x
def by_simpson(f, a=1.0, b=2.0, n=50): # (b-a)/3 Σ[f(a) + 4f(a+h) + 2f(a + h) + 4f(a+2h) + 2f(a+3h) ... + f(b)]
if n % 2 == 1:
n += 1
dx = (b - a) / n
ddx = dx + dx # double dx
x4 = a + dx # 2nd term * 4
x2 = a + ddx # 3rd term * 2
sum = f(a) + f(b) # head + tail
for i in range(1, n - 2, 2) : # exclude head and tail
sum += f(x4) * 4 + f(x2) * 2
x4 += ddx # next
x2 += ddx # next
sum += f(x4) * 4 # last one
return sum * dx / 3
f = lambda x: np.exp(-x**2)
a = 0.0
b = 10.0
n = 15
anser = np.sqrt(np.pi) / 2
jacobi = by_jacobi(f, a, b, n)
newton = by_newton(f, a, b, n)
trapezoid = by_trapezoid(f, a, b, n)
simpson = by_simpson(f, a, b, n)
print(f"∫ :\tjacobi={jacobi} \t, newton={newton} \t, trapezoid={trapezoid}\t, simpson={simpson}\t, compare to {anser}")
print(f"Δ :\t {jacobi - anser} \t, {newton - anser}\t, {trapezoid - anser} \t, {simpson - anser}")