/*
cont1.c 2変数関数の等高線図描画プログラム

(c)2008-2012 和田維作
--------------------------------------------------------------------
●Microsoft Visual Studio .NET 2003、2005、2008 の場合

[ファイル]→[新規作成]→[プロジェクト]を選んで、
[VisualC++]→[Win32]→[Win32プロジェクト]を選び、
プロジェクト名を "<名前を入力してください>" から
"cont1" へ変更する。
その時、V.S 2005、2008 の場合、
[ソリューションのディレクトリを作成]の
チェックボックスをOFFにする。
[OK]を押してプロジェクトを作成する。

Win32アプリケーション ウィザード で左側の
「アプリケーションの設定」タブを開き、
追加のオプションとして、「空のプロジェクト」にチェックを入れて、
「完了」ボタンを押す。

最小化して、新しく作成された cont1 フォルダに cont1.c をコピーする。

最大化して、ソリューションエクスプローラの2行目 cont1 を
クリックして選択する。

[プロジェクト]→[既存項目の追加]で cont1.c をダブルクリックする。
[ビルド]→[構成マネージャ]で[アクティブソリューション構成]を
Releaseにして[閉じる]を押す。

ソリューションエクスプローラの2行目 cont1 を
右クリックしてプロパティ選択する。

[構成プロパティ]の[C/C++]→[コード生成]クリックする。
[ランタイム ライブラリ]を[マルチスレッド(/MT)]に変更する。
[OK]を押す。

続いて[デバッグ]→[デバッグなしで開始]を選んで、[はい]で実行。
--------------------------------------------------------------------
*/
#include <windows.h>
#include <math.h>
#include <process.h>

double f(double x,double y) { return -x*x+y*y*(2.5-y); }

#define INTERVAL 1.0 /* 等高線の間隔 */
#define CENTER_X 0.0 /* ウィンドウの中心を空間上で設定する */
#define CENTER_Y 0.8
#define SCALE    2.0 /* これを半径とした円が入る長方形となる */
#define BASE     0.0 /* 基準の高さ */
#define BACK_COL RGB(255,255,255) /* 背景の色 */
#define ZERO_COL RGB(  0,  0,  0) /* 高さがBASEの等高線の色 */
#define PLUS_COL RGB(  0,255,255) /* BASEより高い等高線の色 */
#define NEGA_COL RGB(  0,  0,255) /* BASEより低い等高線の色 */
#define BGR(A) (((A)&255)<<16|((A)&0xff00)|(A)>>16)

#if defined(MASM32)
 #define _beginthreadex CreateThread
 #define _endthreadex   ExitThread
 #define TH_T DWORD WINAPI
#else
 #define TH_T unsigned WINAPI
#endif

CRITICAL_SECTION Crit;
LONG Tail;

#if _MSC_VER >= 1400
 #define FETCH_INC(i) i=_InterlockedExchangeAdd(&Tail,1);
#else
 #define FETCH_INC(i) \
{ EnterCriticalSection(&Crit); i=Tail++; LeaveCriticalSection(&Crit); }
#endif

__int64**A;
int ASize,*B,BSize,NumOfCpu;
HANDLE*H;
HWND hMain;
HPEN hPen;
HBRUSH hBrs;
int RawMode,BackCol,NegaCol,PlusCol,ZeroCol;
SYSTEM_INFO Info;

volatile double Hy,Hx,Dxy; volatile int Si,Sj,Si2,Sj2;

__int64 IntFunc(double x,double y)
{ return(__int64)floor((f(x,y)-BASE)*(1.0/INTERVAL)); }

TH_T ThFunc(void*arg) {
    int i,j,si2=Si2,sj2=Sj2; __int64*a; double x,y,d=Dxy,hx=Hx,hy=Hy;
    for (i=(int)(size_t)arg;i<si2;) {
        for (a=A[i],y=hy-d*i,x=hx,j=0;j<sj2;j++,x+=d)
        a[j]=IntFunc(x,y);
        FETCH_INC(i);
    }
    _endthreadex(0); return 0;
}

void OnPaint(void) {
    RECT rc; PAINTSTRUCT ps; HDC dc;
    HBITMAP hBitMap,hOld; HDC dcBuffer;
    double x,y,sij; int i,j,k; __int64*d,c;
    dc=BeginPaint(hMain,&ps);
    GetClientRect(hMain,&rc); Sj=rc.right; Si=rc.bottom;
    Sj=((Sj+3)/4)*4;
    Si2=Si+2; Sj2=Sj+2; i=sizeof*A*Si2+sizeof*d*Si2*Sj2;
    if (i>ASize) {
        if (ASize) GlobalFree(A);
        A=(__int64**)GlobalAlloc(GMEM_FIXED,i);
        ASize=i;
    }
    d=(__int64*)(A+Si2); sij=min(Si,Sj); Dxy=2.0*SCALE/(sij-1);
    Hx=CENTER_X-SCALE*Sj/sij-Dxy; Hy=CENTER_Y+SCALE*Si/sij+Dxy;
    for (i=0;i<Si2;i++) { A[i]=d; d+=Sj2; }
    Tail=NumOfCpu;
    if (NumOfCpu>1) {
        for (i=0;i<NumOfCpu;i++)
            H[i]=(HANDLE)_beginthreadex(0,0,ThFunc,(void*)(size_t)i,0,0);
        for (i=0;i<NumOfCpu;i++)
        { WaitForSingleObject(H[i],0x7fffffff); CloseHandle(H[i]); }
    } else {
        for (i=0;i<Si2;i++)
         for (y=Hy-Dxy*i,x=Hx,j=0;j<Sj2;j++,x+=Dxy)
             A[i][j]=IntFunc(x,y);
    }
    i=sizeof*d*Si*Sj;
    if (i>BSize) {
        if (BSize) GlobalFree(B);
        B=(int*)GlobalAlloc(GMEM_FIXED,i);
        BSize=i;
    }
    dcBuffer=CreateCompatibleDC(dc);
    for (i=0;i<Si;i++) for (j=0;j<Sj;j++) {
        c=A[i+1][j+1];
        if (c>A[i+1][j+2]||c>A[i+1][j]||c>A[i+2][j+1]||c>A[i][j+1])
            B[i*Sj+j]=(c<0?NegaCol:c>0?PlusCol:ZeroCol);
        else B[i*Sj+j]=BackCol;
    }
    if (RawMode) hBitMap=CreateBitmap(Sj,Si,1,32,B);
    else         hBitMap=CreateCompatibleBitmap(dc,Sj,Si);
    hOld=(HBITMAP)SelectObject(dcBuffer,hBitMap);
    if (RawMode==0) {
        SelectObject(dcBuffer,hBrs);
        SelectObject(dcBuffer,hPen);
        Rectangle(dcBuffer,0,0,Sj,Si);
        for (k=i=0;i<Si;i++) for (j=0;j<Sj;j++,k++)
         if (B[k]!=BackCol) SetPixel(dcBuffer,j,i,B[k]);
    }
    BitBlt(dc,0,0,Sj,Si,dcBuffer,0,0,SRCCOPY);
    SelectObject(dcBuffer,hOld);
    DeleteObject(hBitMap);
    DeleteDC(dcBuffer);
    EndPaint(hMain,&ps);
}

LRESULT CALLBACK WndProc(HWND hwnd,UINT msg,WPARAM wp,LPARAM lp) {
    switch (msg) {
    case WM_DESTROY: PostQuitMessage(0);          return 0;
    case WM_PAINT:   OnPaint();                   return 0;
    case WM_SIZE:    InvalidateRect(hwnd,0,TRUE); return 0;
    case WM_GETMINMAXINFO:
         ((MINMAXINFO*)lp)->ptMinTrackSize.y=100; return 0;
    default: return DefWindowProc(hwnd,msg,wp,lp);
    }   
}

int WINAPI WinMain(HINSTANCE hinst,HINSTANCE prv,char*arg,int cmd) {
    MSG msg; HDC dc; static WNDCLASSEX wc; static void*dummy;
    dummy=&arg; dummy=&prv; wc.cbSize=sizeof wc;
    GetSystemInfo(&Info); NumOfCpu=Info.dwNumberOfProcessors;
    InitializeCriticalSection(&Crit); 
    if (NumOfCpu>1) H=GlobalAlloc(0,NumOfCpu*sizeof*H);
    wc.lpfnWndProc=WndProc; wc.hInstance=hinst;
    wc.hCursor=LoadCursor(0,IDC_ARROW);
    wc.hbrBackground=GetStockObject(NULL_BRUSH);
    wc.lpszClassName=TEXT("CONT");
    RegisterClassEx(&wc);
    hMain=CreateWindow(TEXT("CONT"),TEXT("CONT"),
        WS_OVERLAPPEDWINDOW,CW_USEDEFAULT,CW_USEDEFAULT,
        CW_USEDEFAULT,CW_USEDEFAULT,0,0,hinst,0);
    dc=GetDC(hMain);
    RawMode=(GetDeviceCaps(dc,BITSPIXEL)==32);
    ReleaseDC(hMain,dc);
    BackCol=(RawMode?BGR(BACK_COL):BACK_COL);
    NegaCol=(RawMode?BGR(NEGA_COL):NEGA_COL);
    PlusCol=(RawMode?BGR(PLUS_COL):PLUS_COL);
    ZeroCol=(RawMode?BGR(ZERO_COL):ZERO_COL);
    hBrs=CreateSolidBrush(BackCol);
    hPen=CreatePen(1,PS_SOLID,BackCol);
    ShowWindow(hMain,cmd);
    while (GetMessage(&msg,0,0,0)) DispatchMessage(&msg);
    if (NumOfCpu>1) GlobalFree(H);
    DeleteCriticalSection(&Crit);
    DeleteObject(hBrs); DeleteObject(hPen);
    GlobalFree(A); return(int)msg.wParam;
}