solveview.cpp

经典四阶龙哥库塔法解方程组

// SolveView.cpp : implementation of the CSolveView class
//
#include <math.h>
#include "stdafx.h"
#include "Solve.h"
#include <iostream>
#include "SolveDoc.h"
#include "SolveView.h"

#ifdef _DEBUG
#define new DEBUG_NEW
#undef THIS_FILE
static char THIS_FILE[] = __FILE__;
#endif

/////////////////////////////////////////////////////////////////////////////
// CSolveView

IMPLEMENT_DYNCREATE(CSolveView, CView)

BEGIN_MESSAGE_MAP(CSolveView, CView)
	//{{AFX_MSG_MAP(CSolveView)
	ON_COMMAND(IDM_SOLVE, OnSolve)
	ON_COMMAND(ID_OUTPUT, OnOutput)
	//}}AFX_MSG_MAP
	// Standard printing commands
	ON_COMMAND(ID_FILE_PRINT, CView::OnFilePrint)
	ON_COMMAND(ID_FILE_PRINT_DIRECT, CView::OnFilePrint)
	ON_COMMAND(ID_FILE_PRINT_PREVIEW, CView::OnFilePrintPreview)
END_MESSAGE_MAP()

/////////////////////////////////////////////////////////////////////////////
// CSolveView construction/destruction

CSolveView::CSolveView()
{
	// TODO: add construction code here
	bSolved = FALSE;
	penGreen.CreatePen(PS_SOLID,1,RGB(255,255,0));
	penRed.CreatePen(PS_SOLID,1,RGB(255,0,0));
	penBlue.CreatePen(PS_SOLID,1,RGB(0,0,255));
	penBlack.CreatePen(PS_SOLID,1,RGB(0,0,0));
}

CSolveView::~CSolveView()
{
	DeleteObject(&penRed);
	DeleteObject(&penGreen);
	DeleteObject(&penBlue);
	DeleteObject(&penBlack);
}

BOOL CSolveView::PreCreateWindow(CREATESTRUCT& cs)
{
	// TODO: Modify the Window class or styles here by modifying
	//  the CREATESTRUCT cs

	return CView::PreCreateWindow(cs);
}

/////////////////////////////////////////////////////////////////////////////
// CSolveView drawing

void CSolveView::OnDraw(CDC* pDC)
{
	CSolveDoc* pDoc = GetDocument();
	ASSERT_VALID(pDoc);
	// TODO: add draw code for native data here
	CRect rect;
   	GetClientRect(&rect);

	if(bSolved)
		DrawSolves(&m_memDC);
	/////////////////////////////////
	pDC->BitBlt(0,0,rect.Width(),rect.Height(),&m_memDC,0,0,SRCCOPY);
}

/////////////////////////////////////////////////////////////////////////////
// CSolveView printing

BOOL CSolveView::OnPreparePrinting(CPrintInfo* pInfo)
{
	// default preparation
	return DoPreparePrinting(pInfo);
}

void CSolveView::OnBeginPrinting(CDC* /*pDC*/, CPrintInfo* /*pInfo*/)
{
	// TODO: add extra initialization before printing
}

void CSolveView::OnEndPrinting(CDC* /*pDC*/, CPrintInfo* /*pInfo*/)
{
	// TODO: add cleanup after printing
}

/////////////////////////////////////////////////////////////////////////////
// CSolveView diagnostics

#ifdef _DEBUG
void CSolveView::AssertValid() const
{
	CView::AssertValid();
}

void CSolveView::Dump(CDumpContext& dc) const
{
	CView::Dump(dc);
}

CSolveDoc* CSolveView::GetDocument() // non-debug version is inline
{
	ASSERT(m_pDocument->IsKindOf(RUNTIME_CLASS(CSolveDoc)));
	return (CSolveDoc*)m_pDocument;
}
#endif //_DEBUG

/////////////////////////////////////////////////////////////////////////////
// CSolveView message handlers
void CSolveView:: InitMemDC()
{
//	CDC *pDC = GetDC();
	CRect rect;
	GetClientRect(&rect);
	CClientDC dc(this);
	m_memDC.CreateCompatibleDC(&dc);

	CBitmap bitmap;
	bitmap.CreateCompatibleBitmap(&dc,rect.Width(),rect.Height());

	m_memDC.SelectObject(&bitmap);
	m_memDC.FillSolidRect(&rect,RGB(168,168,168));
//	m_memDC.BitBlt(0,0,rect.Width(),rect.Height(),&dc,0,0,SRCCOPY);
}

void CSolveView::OnInitialUpdate() 
{
	CView::OnInitialUpdate();
	
	// TODO: Add your specialized code here and/or call the base class

}

void CSolveView::DrawSolves(CDC *pDC)
{
//	pDC->SetMapMode(MM_LOMETRIC);
	int i;
	CString s;
	pDC->SetTextColor(RGB(0,0,0));
	for(i = 0; i<n ; i++)
	{
		s.Format("x1[%f]=%f,  x2[%f]=%f,  x3[%f]=%f",t[i],x1[i],t[i],x2[i],t[i],x3[i]);
		pDC->TextOut(0+20,0+i*20,s);
	}

/////////////////////////////////////////////////////
	int xoff = 600;
	int yoff = 600;

	oldPen = pDC->SelectObject(&penRed);
	pDC->MoveTo(750,20);
	pDC->LineTo(900,20);
	pDC->TextOut(730,10," x1 ");
	oldPen = pDC->SelectObject(&penGreen);
	pDC->MoveTo(750,40);
	pDC->LineTo(900,40);
	pDC->TextOut(730,30," x2 ");
	oldPen = pDC->SelectObject(&penBlue);
	pDC->MoveTo(750,60);
	pDC->LineTo(900,60);
	pDC->TextOut(730,50," x3 ");

//////////////////////////////////////////////////////////////////
	oldPen = pDC->SelectObject(&penBlack);
	pDC->MoveTo(0+xoff,0+yoff);
	pDC->TextOut(0+xoff-10,0+yoff-10,"0");

	pDC->LineTo(0+xoff+380,0+yoff);
	pDC->LineTo(0+xoff+380-5,0+yoff-5);
	pDC->MoveTo(0+xoff+380,0+yoff);
	pDC->LineTo(0+xoff+380-5,0+yoff+5);
	pDC->TextOut(0+xoff+380-5,0+yoff+5," + t ");

	pDC->MoveTo(0+xoff,0+yoff);
	pDC->LineTo(0+xoff,0+yoff-550);
	pDC->LineTo(0+xoff-5,0+yoff-550+5);
	pDC->MoveTo(0+xoff,0+yoff-550);
	pDC->LineTo(0+xoff+5,0+yoff-550+5);
	pDC->TextOut(0+xoff+5,0+yoff-550-5," + x ");

/*	pDC->MoveTo(0+xoff,0+yoff);
	pDC->LineTo(0+xoff,0+yoff+300);
	pDC->LineTo(0+xoff-5,0+yoff+300-5);
	pDC->MoveTo(0+xoff,0+yoff+300);
	pDC->LineTo(0+xoff+5,0+yoff+300-5);
	pDC->TextOut(0+xoff+5,0+yoff+300-5," - x ");
*/


//--将坐标绕x轴反转了，这才是正确的，也可以通过设置映射模式改变--------------------
/*	pDC->SelectObject(&penRed);
	pDC->MoveTo(t[0]+xoff,-x1[0]+yoff);
	for(i=0;i<n;i++)
	{
		pDC->LineTo(t[i]+xoff,-x1[i]+yoff);
	}
	pDC->SetTextColor(RGB(255,0,0));
	pDC->TextOut(t[n]+xoff,-x1[n]+yoff," x1(t) ");


	pDC->SelectObject(&penGreen);
	pDC->MoveTo(t[0]+xoff,-x2[0]+yoff);
	for(i=0;i<n;i++)
	{
		pDC->LineTo(t[i]+xoff,-x2[i]+yoff);
	}
	pDC->SetTextColor(RGB(0,255,0));
	pDC->TextOut(t[n]+xoff,-x2[n]+yoff," x2(t) ");


	pDC->SelectObject(&penBlue);
	pDC->MoveTo(t[0]+xoff,-x3[0]+yoff);
	for(i=0;i<n;i++)
	{
		pDC->LineTo(t[i]+xoff,-x3[i]+yoff);
	}
	pDC->SetTextColor(RGB(0,0,255));
	pDC->TextOut(t[n]+xoff,-x3[n]+yoff," x3(t) ");
*/


	oldPen = pDC->SelectObject(&penRed);
	pDC->MoveTo(10*t[0]+xoff,-x1[0]*0.7+yoff);
	for(i=0;i<n;i++)
	{
		pDC->LineTo(10*t[i]+xoff,-x1[i]*0.7+yoff);
	}
	pDC->SetTextColor(RGB(255,0,0));
	pDC->TextOut(10*t[n]+xoff,-x1[n]*0.7+yoff," x1(t) ");


	oldPen = pDC->SelectObject(&penGreen);
	pDC->MoveTo(10*t[0]+xoff,-x2[0]*0.7+yoff);
	for(i=0;i<n;i++)
	{
		pDC->LineTo(10*t[i]+xoff,-x2[i]*0.7+yoff);
	}
	pDC->SetTextColor(RGB(0,255,0));
	pDC->TextOut(10*t[n]+xoff,-x2[n]*0.7+yoff," x2(t) ");


	oldPen = pDC->SelectObject(&penBlue);
	pDC->MoveTo(10*t[0]+xoff,-x3[0]*0.7+yoff);
	for(i=0;i<n;i++)
	{
		pDC->LineTo(10*t[i]+xoff,-x3[i]*0.7+yoff);
	}
	pDC->SetTextColor(RGB(0,0,255));
	pDC->TextOut(10*t[n]+xoff,-x3[n]*0.7+yoff," x3(t) ");
	pDC->SelectObject(oldPen);

//-------------------------------------------------------------------------
}


///////////////////////////////////////////////////////////////////////////////
float  CSolveView::f1(float t,float x1,float x2,float x3)	//x1`=f1(t,x1,x2,x3)
{
	float r;
//	return r= -0.008*x1-0.1*x2;
	return r = -0.001*x1*x2;
}
float  CSolveView::f2(float t,float x1,float x2,float x3)	//x2`=f2(t,x1,x2,x3)
{
	float r;
//	return r= -0.004*t;
	return r = 0.001*x1*x2 - 0.072*x2;
}
float  CSolveView::f3(float t,float x1,float x2,float x3)	//x3`=f3(t,x1,x2,x3)
{
	float r;
//	return r= -0.005*x2+0.2*x1;
	return r = 0.072*x3;
}

void CSolveView::RungeKutta(float t0,float x10, float x20, float x30, float step, int num)
{
	t[0] = t0;
	x1[0] = x10;
	x2[0] = x20;
	x3[0] = x30;
	h = step;
	n = num;

	float k1,k2,k3,k4;
	float l1,l2,l3,l4;
	float m1,m2,m3,m4;

	for(int i=0; i<n; i++ )
	{
		k1 = h*f1( t[i],x1[i],x2[i],x3[i] );
		l1 = h*f2( t[i],x1[i],x2[i],x3[i] );
		m1 = h*f3( t[i],x1[i],x2[i],x3[i] );

		k2 = h*f1( t[i]+h/2,x1[i]+k1/2,x2[i]+l1/2,x3[i]+m1/2 );
		l2 = h*f2( t[i]+h/2,x1[i]+k1/2,x2[i]+l1/2,x3[i]+m1/2 );
		m2 = h*f3( t[i]+h/2,x1[i]+k1/2,x2[i]+l1/2,x3[i]+m1/2 );

		k3 = h*f1( t[i]+h/2,x1[i]+k2/2,x2[i]+l2/2,x3[i]+m2/2 );
		l3 = h*f2( t[i]+h/2,x1[i]+k2/2,x2[i]+l2/2,x3[i]+m2/2 );
		m3 = h*f3( t[i]+h/2,x1[i]+k2/2,x2[i]+l2/2,x3[i]+m2/2 );
		
		k4 = h*f1( t[i]+h,x1[i]+k3,x2[i]+l3,x3[i]+m3 );
		l4 = h*f2( t[i]+h,x1[i]+k3,x2[i]+l3,x3[i]+m3 );
		m4 = h*f3( t[i]+h,x1[i]+k3,x2[i]+l3,x3[i]+m3 );

		x1[i+1] = x1[i] + (k1+2*k2+2*k3+k4)/6;
		x2[i+1] = x2[i] + (l1+2*l2+2*l3+l4)/6;
		x3[i+1] = x3[i] + (m1+2*m2+2*m3+m4)/6;

		t[i+1] = t[i] + h;
	}

	bSolved = TRUE;
}


void CSolveView::OnSolve() 
{
	// TODO: Add your command handler code here
	InitMemDC();
	RungeKutta(0, 720, 10, 70, 0.1, 1999);
	Invalidate();
}

void CSolveView::OnOutput() 
{
	// TODO: Add your command handler code here

	if(!bSolved)
		return;

	CFileDialog fileDlg(FALSE);
	fileDlg.m_ofn.lpstrTitle="数据保存";
	fileDlg.m_ofn.lpstrFilter="Text Files(*.txt)\0*.txt\0All Files(*.*)\0*.*\0\0";
	
	if(IDOK==fileDlg.DoModal())
	{
		CString s;
		CFile file(fileDlg.GetFileName(),CFile::modeCreate | CFile::modeWrite);//根据文件名构造
	//	file.Write("经典四阶龙哥库塔法解微分方程组\r\n",strlen("经典四阶龙哥库塔法解微分方程组"));
		file.Write("t           x(1)           x(2)          x(3) \r\n",strlen("t           x(1)           x(2)          x(3)"));

		for(int i=0;i<2000;i++)
		{
			s.Format("%f,  %f,  %f,  %f, \r\n",t[i],x1[i],x2[i],x3[i]);
			file.Write(s,strlen(s));
		}
		file.Close();
	}

}