SakNgoi

The dynamics of Plasmodium falciparum during the expansion phase of the asexual stage of infection

https://www.wolframcloud.com/objects/sompob/models/NJWandImmune.nb

https://pantip.com/topic/38221485/

code สำหรับปัญหานี้

https://www.wolframcloud.com/objects/sompob/test/pantip_01Nov18.nb

Manipulate[
Show[
  Graphics[{EdgeForm[Pink], FaceForm[Blue], 
    Rectangle[{1 - Cos[-15 Degree], 0}, {x, 
      x - (1 - Cos[-15 Degree])}]}],
  Graphics[{EdgeForm[LightBlue], FaceForm[Black], 
    Rotate[Rectangle[{0, 0}], -15 Degree, {1, 0}],
    PlotRange -> {{-2, 2}, {-2, 2}} }],
  Graphics[
   Text["Area of the Big Square:" <> 
     ToString[(x - (1. - Cos[15 Degree]))^2/(0.5*Sin[15 Degree]*
         Cos[15 Degree])], {0.25, 1.5}]]
   ],
 {{x, 1.259, "side"}, 0.1, 2, 0.001}
 ]

ทดลองเล่น WolframCloud

เจ๋งมากกกกกก เขียนสั้น deploy ก็ง่าย (เชียร์ให้คนมาใช้กันเยอะๆ 555)

ลองดูตัวอย่างที่ผมเขียนนี้ครับ  ไอเดียก็คือว่าให้มันเขียนภาพเป็นไฟล์ pdf ด้วยคำสั่ง ggsave แล้วเรียกกลับมาใน Mathematica ด้วย Import  โค้ดที่เขียนนี้ดัดแปลงมาจาก https://mathematica.stackexchange.com/questions/16726/how-to-get-a-plot-generated-by-r-returned-in-an-output-cell

โค้ดที่เขียนด้วยMathematicaหรือภาษา Wolfram นั้นมันจะมีวิธีที่ช่วยให้มันรันได้เร็วขึ้นอยู่หลายวิธีครับ ขึ้นกับปัญหาและวิธีการเขียน ที่จะโชว์ให้ดูนี้ผมก็ใช้เทคนิคง่ายๆด้วยการทำ parallel ด้วยคำสั่ง ParallelTable และ compile ด้วยคำสั่ง Compile ครับ  โค้ดตัวอย่างที่จะเอามาลองนี้เป็นโค้ดที่ใช้วาดรูป Mandelbrot set ครับ

แบบที่ 1 เป็นโค้ดเริ่มต้นที่ยังไม่มีการทำให้มันเร็วครับเขียนโดยใช้ loop จากคำสั่ง Table

AbsoluteTiming[
 test1 = Block[{i, x, p}, 
    Table[i = 0; x = 0. I; p = r + I c; 
     While[Abs@x <= Sqrt[2] && i < 9^3, x = x^2 + p; ++i]; 
     Tanh[Power[i/9^3, (7)^-1]], {c, -1, 1, .01}, {r, -2, 1, .01}]];]

(* เวลาที่ใช้คือ 77.4437 วินาที *)

แบบที่ 2 ในแบบนี้ผมทดลองให้มันรันแบบขนานโดยใช้คำสั่ง ParallelTable เครื่องที่ผมใช้นี้มันมี 4 cores ครับ เวลาก็ลดลง ถ้าคิดง่ายๆก็ใช้เวลาเหลือเกือบๆ 1 ใน 4

AbsoluteTiming[
 test2 = Block[{i, x, p}, 
    ParallelTable[i = 0; x = 0. I; p = r + I c; 
     While[Abs@x <= Sqrt[2] && i < 9^3, x = x^2 + p; ++i]; 
     Tanh[Power[i/9^3, (7)^-1]], {c, -1, 1, .01}, {r, -2, 1, .01}]];]

(* เวลาที่ใช้คือ 23.8938 วินาที *)

แบบที่ 3  แบบนี้ก็คือเอาแบบที่ 1 มา compile เป็น binary code ด้วยคำสั่ง Compile เวลาก็เร็วขึ้นอย่างเห็นได้ชัด

test3 = Compile[{}, 
  Block[{i, x, p}, 
   Table[i = 0; x = 0. I; p = r + I c; 
    While[Abs@x <= Sqrt[2] && i < 9^3, x = x^2 + p; ++i]; 
    Tanh[Power[i/9^3, (7)^-1]], {c, -1, 1, .01}, {r, -2, 1, .01}]]]

AbsoluteTiming[test3[];]

(* เวลาที่ใช้คือ 2.68316 วินาที *)

แบบที่ 4  แบบนี้ก็คือแบบที่ 3 ที่เอามารันแบบขนานด้วยคำสั่ง ParallelTable ซึ่งเวลาที่ใช้ก็เร็วขึ้นอีกกว่าแบบที่ 3

test4 = Compile[{{c, _Real}, {r, _Real}},
   Module[{i, x, p},
    i = 0; x = 0. I; p = r + I c; 
    While[Abs@x <= Sqrt[2] && i < 9^3, x = x^2 + p; ++i]; 
    Tanh[Power[i/9^3, (7)^-1]]
    ]
   ];

ParallelTable[
   test4[c, r], {c, -1, 1, .01}, {r, -2, 1, .01}]; // AbsoluteTiming

(* เวลาที่ใช้คือ 1.04455 วินาที *)

แบบที่ 5 แบบนี้เหมือนกับแบบที่ 4 แต่แทนที่จะให้มัน compile เป็น binary code โดย Mathematica เอง ก็ให้มันแปลงเป็นภาษา C ด้วย option CompilationTarget -> “C” เลยโดยใช้ compiler อย่าง gcc ครับ เวลาก็เร็วขึ้นอีกประมาณ 60% จากแบบที่ 4 ครับ

Needs["CCompilerDriver`GenericCCompiler`"]

$CCompiler = {"Compiler" -> GenericCCompiler, 
   "CompilerInstallation" -> "C:/Rtools/mingw_64", 
   "CompilerName" -> "x86_64-w64-mingw32-gcc.exe"};


test5 = Compile[{{c, _Real}, {r, _Real}},
   Module[{i, x, p},
    i = 0; x = 0. I; p = r + I c; 
    While[Abs@x <= Sqrt[2] && i < 9^3, x = x^2 + p; ++i]; 
    Tanh[Power[i/9^3, (7)^-1]]
    ], CompilationTarget -> "C"];

ParallelTable[
   test5[c, r], {c, -1, 1, .01}, {r, -2, 1, .01}]; // AbsoluteTiming

(* เวลาที่ใช้ 0.393584 วินาที *)

แนะนำ Link นี้ครับสำหรับผู้สนใจใช้งาน Mathematica ขุมทรัพย์อีกแห่งเลยครับ 🙂

https://stackoverflow.com/questions/4198961/what-is-in-your-mathematica-tool-bag

เทคนิคเจ๋งๆเพียบเลยครับ

ใน Rtools จะมี compiler ของ gcc (mingw) มาด้วยแล้วทั้งที่เป็นแบบ 32 bits และ 64 bits หากอยากจะเอาไปใช้ใน Mathematica ก็ต้องเรียกผ่าน CCompilerDriver`GenericCCompiler` ครับ แล้วเพียงเซ็ต Path ของ gcc จาก Rtools นี้ให้ถูก เช่น

Needs["CCompilerDriver`GenericCCompiler`"]

ทดลองเรียกใช้งาน

greeter = CreateExecutable[StringJoin[
 "#include <stdio.h>\n",
 "int main(){\n",
 " printf(\"Hello MinGW-w64 world.\\n\");\n",
 "}\n"],
 "helloworld", "Compiler" -> GenericCCompiler, 
 "CompilerInstallation" -> "C:/Rtools/mingw_64", 
 "CompilerName" -> "x86_64-w64-mingw32-gcc.exe"]

Import["!\""<>greeter<>"\"","Text"]
Hello MinGW-w64 world.

ที่นี้ถ้าอยากจะให้ Mathematica มันเรียกใช้เจ้า mingw 64 จาก Rtools นี้ตลอด อย่างเช่นคำสั่ง Compile ก็สามารถทำได้โดยเซ็ตค่าที่มันเกี่ยวข้องอย่างเช่น Path กับตัวแปร $CCompiler ได้เลยครับ

$CCompiler = {"Compiler" -> GenericCCompiler, 
 "CompilerInstallation" -> "C:/Rtools/mingw_64", 
 "CompilerName" -> "x86_64-w64-mingw32-gcc.exe"};

f = Compile[{x, y}, Sqrt[x^2 + y^2], CompilationTarget -> "C"]

Table[{x, f[x, 5/4 x]}, {x, 0, 6, 0.5}]

ปกติผมจะใช้ Mathematica เป็นหลักในการทำงานโดยมากจะเป็นพวกแบบจำลองคณิตศาสตร์ที่ตองการเห็นว่าผลลัพธ์จะเปลี่ยนแปลงไปอย่างไรถ้าตัวแปรที่สนใจบางตัวมีการเปลี่ยนแปลง เพิ่มขึ้นหรือลดลง ซึ่งแน่นอนผมใช้พวกคำสั่ง Manipulate หรือ Dynamic เป็นหลัก แต่มาระยะหลังนี้เริ่มมีการใช้ R ในที่ทำงานกันเยอะ ผมเลยเขียนโมเดลของงานที่ทำใน R ไปพอสมควร แต่มันก็ติดปัญหาว่ามันช้ามากกับหลายโมเดลที่ทำอยู่ถึงแม้จะหันไปใช้ Rcpp แต่ก็รู้สึกว่าเสียเวลาเขียน wrapper functions ห่อพวกมันอีกทีใน R แถมเวลาส่งงานให้ Prof ทั้งหลายก็ติดเรื่อง compiler อีก แถมต้องมาเขียน shiny อีกดูวุ่นวายไปหมด สุดท้ายกลับมาตายรังที่ Mathematica เหมือนเดิม เขียนโมเดลเสร็จ อยากจะส่งให้ชาวบ้านดูก็ใช้เพียงแค่คำสั่ง Deploy ไปบนcloud เท่านั้นจบเลย …แต่มันก็แลกมาด้วยค่าใช้จ่ายพอสมควร

ผมเลยหาทางเลือกอื่นๆที่มันประหยัดกับบางงานกับโมเดลที่ไม่ดูซับซ้อนและใช้เวลาทำไม่นาน หันซ้ายหันขวาก็มาจบที่ C# หลังจากดูแล้วก็มี library ที่พอช่วยงานได้อย่าง math.net แถม UI ก็ทำง่ายเพียงลากแปะเท่านั้น เลยทดลองใช้ดู ผลปรากฏว่าเป็นที่น่าพอใจ  ถึงจะดูยุ่งยากในส่วน UI แต่ก็รู้สึกว่าไม่ได้เลวร้ายอะไร ส่งงานก็ง่ายแถมก็ไม่ต้องติดตั้งอะไรให้ดูวุ่นวาย

ตัวอย่างงงานที่ทำที่พอให้ดูได้ครับ

สนใจอยากทดลองเล่นก็ download ได้เลยครับ

พอดีว่าลงWindows 7 64bit กับ Mathematica 11.1.1 ใหม่ แล้วเจอปัญหา  ตามภาพด้านล่างตอนใช้ CUDALink พอเช็คดูเวอร์ชั่นของ nvidia driver ที่ใช้กับการ์ดจอ (Quadro 2000)ด้วยคำสั่ง CUDADriverVersion[] พบว่าเป็นรุ่นเก่ามากคือ27x.x  เลยทำการ update เป็นตัวใหม่กว่าคือ 377.55จากนั้นก็ลองใช้ QUDAQ[] ใหม่เช็คดูปรากฎว่าไม่มีปัญหาแล้วสามาถใช้งานได้ตามปกติ

เห็นมีคนพูดถึงปัญหานี้ในกลุ่มผู้ใช้ R ของไทยเราเลยเอามาเขียนโปรแกรมทำบ้างแต่ใช้ภาษา Wolfram แทนเพราะจะเขียนได้สั้นมาก

โปรแกรมที่เขียนนี้จะวาดเส้นแนวนอนที่มีระยะห่างตามที่กำหนดและก็สุ่มวางเข็มตามขนาดและจำนวนที่ต้องการลงไปแล้วนับดูว่าที่กี่เล่มที่ทับเส้น จากนั้นก็คำนวณประมาณค่า pi จากสูตรที่เขียนไว้ที่ https://en.wikipedia.org/wiki/Buffon%27s_needle ซึ่งก็คือ

โดย l คือความยาวเข็ม, n  คือจำนวนเข็ม, t คือระยะห่างระหว่างเส้นบรรทัด และ h คือจำนวนเข็มที่ทับเส้นบรรทัด

สำหรับวาดเส้นของเข็มและเส้นบรรทัด

genNeedle[{x_, y_}, θ_, len_] :=
Line[{{(x + 0.5 len Cos[θ]), (y + 0.5 Sin[θ])}, {(x –
0.5 len Cos[θ]), (y – 0.5 len Sin[θ])}}]

horizontalLine[diff_, len_, {x0_, y0_}, n_] :=
Table[Line[{{x0, y0 + (i diff)}, {len, y0 + (i diff)}}], {i, n}]

สำหรับประมาณค่า pi

estimatedPi[needlelength_, nneedles_, linespace_, ncross_] :=
2.*needlelength*nneedles /(linespace*ncross)

สำหรับ simulation

run[nl_, nfl_, ll_, diff_] :=
Module[{fls, ls, outls, pnts, ncross},
fls = horizontalLine[diff, 10, {0, 0}, nfl];
ls = Table[
genNeedle[{RandomReal[{0, 10}], RandomReal[{0, nfl*diff}]},
RandomReal[π], ll], {nl}];
outls =
Table[Solve[{x, y} ∈ RegionIntersection[fls[[i]], #], {x,
y}] & /@ ls, {i, Length@fls}];
pnts = Point[{x, y}] /. Flatten[outls, 2];
ncross = Length@Flatten[outls, 2];
Graphics[{fls, ls, Red, PointSize[Medium],
If[Length@Flatten[outls, 2] > 0, pnts]},
PlotLabel ->
Style["จำนวนเข็มที่ทับเส้น:" <> ToString@ncross <>
" π≈" <>
ToString@estimatedPi[ll, nl, diff, ncross], FontSize -> 20,
Bold]]];

ทดลองเล่นดู โดยใช้เข็ม 100 เล่ม และกำหนดให้มีเส้นบรรทัด 10 เส้น โดยระยะห่างระหว่างบรรทัดกับความยาวเข็มเท่ากัน

nlines = 10;
nneedles = 100;
linelength = 1;
linespace = 1;
run[nneedles, nlines, linelength, linespace]

สำหรับmathematica notebook สามารถdownloadได้จากhttps://github.com/slphyx/Buffon-s-Needle