e^t

哇呀呀……

工作是不是因为无聊？学习是不是因为无聊？介个不晓得里～～

话说卖萌可以买到妹纸吗？In urgent need ing~

有图有代码

从一个如下格式的文件中读取球的坐标和半径，然后画出来，使用了vtk库

n

x1 y1 z1 r1

x2 y2 z2 r2

…

xn yn zn rn

编译命令如下：

gcc example.cpp -I/usr/include/vtk-5.4 -L/usr/lib -lvtkCommon -lvtkViews -o example

#include <vtkSphereSource.h>
#include <vtkSmartPointer.h>
#include <vtkAppendPolyData.h>
#include <vtkCamera.h>
#include <vtkPolyDataMapper.h>
#include <vtkActor.h>
#include <vtkRenderer.h>
#include <vtkRenderWindow.h>
#include <vtkProperty.h>
#include <vtkRenderWindowInteractor.h>
#include <vtkTransform.h>
#include <vtkDepthSortPolyData.h>

vtkSmartPointer<vtkAppendPolyData> GenerateSpheres(int ntheta,int nphi)
{
    vtkSmartPointer<vtkAppendPolyData> appendData =
        vtkSmartPointer<vtkAppendPolyData>::New();

    ifstream sphereFile("spheres.3D");
    int nSphere;
    sphereFile >> nSphere;
    std::cout << "Number of spheres: " << nSphere << std::endl;
    double coordSphere[nSphere][4];
    for(int i = 0;i<nSphere;i++)
        for(int j = 0;j<4;j++)
            sphereFile >> coordSphere[i][j];

    for(int i = 0;i<nSphere;i++)
    {
        vtkSmartPointer<vtkSphereSource> sphereSource =
            vtkSmartPointer<vtkSphereSource>::New();
        sphereSource->SetThetaResolution(ntheta);
        sphereSource->SetPhiResolution(nphi);
        sphereSource->SetRadius(coordSphere[i][3]);

        sphereSource->SetCenter(coordSphere[i][0],
                coordSphere[i][1],
                coordSphere[i][2]);
        sphereSource->Update();
        appendData->AddInputConnection(sphereSource->GetOutputPort());
    }
    return appendData;
}

int main (int argc, char *argv[])
{
  int theta = 16;
  int phi = 16;

  // Generate a translucent sphere poly data set that partially overlaps:
  vtkSmartPointer<vtkAppendPolyData> translucentGeometry =
    GenerateSpheres(theta, phi);

  // generate a basic Mapper and Actor
  vtkSmartPointer<vtkPolyDataMapper> mapper =
    vtkSmartPointer<vtkPolyDataMapper>::New();
  mapper->SetInputConnection(translucentGeometry->GetOutputPort());

  vtkSmartPointer<vtkActor> actor =
    vtkSmartPointer<vtkActor>::New();
  actor->SetMapper(mapper);
  actor->GetProperty()->SetOpacity(1); //change to make translucent !!!
  actor->GetProperty()->SetColor(1, 0, 0);
  actor->RotateX(0); // put the objects in a position where it is easy to see
                       // different overlapping regions

  // Create the RenderWindow, Renderer and RenderWindowInteractor
  vtkSmartPointer<vtkRenderer> renderer =
    vtkSmartPointer<vtkRenderer>::New();
  vtkSmartPointer<vtkRenderWindow> renderWindow =
    vtkSmartPointer<vtkRenderWindow>::New();
  renderWindow->AddRenderer(renderer);
  vtkSmartPointer<vtkRenderWindowInteractor> renderWindowInteractor =
    vtkSmartPointer<vtkRenderWindowInteractor>::New();
  renderWindowInteractor->SetRenderWindow(renderWindow);

  // Add the actors to the renderer, set the background and size
  renderer->AddActor(actor);
  renderer->SetBackground(1, 1, 1);
  renderWindow->SetSize(600, 400);

  // Setup view geometry
  renderer->ResetCamera();
  renderer->GetActiveCamera()->Zoom(2.2); // so the object is larger

  int success = EXIT_SUCCESS;

  // Initialize interaction
  renderWindowInteractor->Initialize();

  // Start interaction
  renderWindowInteractor->Start();

  return success;
}

 

vtk的一些基础特性（http://apps.hi.baidu.com/share/detail/23071009）：

VTk是在三维函数库OpenGL的基础上，采用面向对象的设计方法发展起来的。它有两种不同的方式：图像模型和可视化模型。

1、图像模型：支持3D几何数据绘制，3D体数据绘制，2D几何文字，图像绘制。

采用了3D图形系统简单易用的特点，同时也采用了图形用户接口的方法。整个图形模型表现3D图形系统的本质特点，主要有9类基本对象：渲染控制器，渲染窗口，渲染器，灯光，照相机，角色，特性，映射，变换。

1）渲染控制器定义与设备无关的坐标计算方法，并创建渲染窗口。

2）渲染窗口管理显示设备上的窗口，一个或多个绘制方法可在渲染窗口上创建一个场景，渲染窗口是用户图形界面，其中包括设置渲染窗口的大小，产生立体显示效果等方法。

3）渲染器是管理光源照相机和绘制对象等的位置、属性等，提供了世界坐标系，观察坐标系和现实坐标系之间的转换。

4）灯光可在场景中照亮绘制对象，可通过调用参数改变控制灯光的状态、照射角度，照射强度，颜色等，并支持点光源和平行光源。

5）照相机是定义观察者的位置，聚焦点和其他相关属性，参数可有调用者根据需要设置。

6）（相当于OpenGL的模型变换）角色代表渲染场景中的绘制对象实体，通过参数调节可以设置角色的位置，方向，渲染特性，引用，纹理映射等属性，并可对角色进行缩放。

7）属性是说明几何物体的一些特性，实现三维图形真实感。

8）映射指定了渲染数据和图形库中基本图元之间的关系，一个或多个角色可以使用相同的映射，有多个参数可对其进行控制。

9）变换是一个放置4×4变换矩阵的堆栈，可以进行各种操作。

2、可视化模型：利用数据流程模型，由过程对象和数据对象组成。过程对象包括可视化流程的模块及算法，数据对象包括数据处理及当数据在网络中流动时对数据进行的操作。过程对象又分为：源对象，过滤器，映射。源对象是可视化流程的起点，源对象包括从文件中读入及程序内部产生的数据，过滤器接收从源对象来的输入数据，处理数据及输出数据，映射指定了基本图元与数据之间的接口，接收过滤输出的数据，并把数据映射为基本图元。

 

References:

[1] http://www.vtk.org/Wiki/VTK/Examples/Cxx/Visualization/CorrectlyRenderingTranslucentGeometry

[2] The official vtk documentation

没事儿熟悉一下C语言，林火模型其实物理学家用来研究SOC(self-organized criticality)，也就是自组织临界性的玩意儿。说到其使用价值，把一个森林简化到这种地步，值得怀疑啊。下面用了狗血的周期边界条件（敢问编过科学计算程序的人，谁没有啃过面包圈？！），初始条件是全部是森林（100*100的网格）。图一表示几个闪电打下了之后，森林鸭梨很大，纷纷dead；图二是图一中的强大火势快没了，因为森林已经七零八落了；图三说明森林恢复力惊人，继续和大火做起了永不停息的斗争……

下面是C语言代码，直接gcc就行了。上面的图我用Mathematica画的。差点忘了，秀一下森林人口变化（刚开始是一万）：

顺便表示对wordpress排版无力吐槽……

#include <stdio.h>

#include <stdlib.h>

#define p_fire 0.001

#define p_grow 0.02

#define ISFIRE 10

#define ISNULL 0

#define ISTREE 1

double ran1()

{

    return 1.0*rand()/RAND_MAX;

}

int main(int argc, const char *argv[])

{

    const int NX = 100;

    const int NY = 100;

    int forest[NX][NY];

    int forest_n[NX][NY];

    int i,j;

    int ide,idw,jdn,jds,sum;

    int n;

    FILE *pFile;

    pFile = fopen("forestfire.dat","w");

    // initialization

    for(i=0;i

        for(j=0;j

            forest[i][j]=1;

            //if(ran1() < 0.818) forest[i][j] = 0;

            //else forest[i][j] = 1;

    // save the initial forest

    for(i=0;i

        for(j=0;j

            fprintf(pFile,"%2d",forest[i][j]);

        fprintf(pFile,"\n");

    }

    for(n=1;n

    {

        for(i=0;i

        {

            for(j=0;j

            {

                ide = i + 1 > NX - 1 ? 0 : i + 1;

                idw = i - 1 < 0 ? NX - 1 : i - 1;

                jdn = j + 1 > NY - 1 ? 0 : j + 1;

                jds = j - 1 < 0 ? NY - 1 : j - 1;

                sum = forest[i][jdn] + forest[i][jds] + forest[ide][j] + forest[idw][j] +

                    forest[ide][jdn] + forest[ide][jds] + forest[idw][jdn] + forest[idw][jds];

                if(forest[i][j] == ISTREE)

                    if(sum>10) // neighbour is on fire!

                        forest_n[i][j] = ISFIRE;

                    else // no neighbour is on fire, lightning sets in however

                        forest_n[i][j] = (ran1()                else if(forest[i][j] == ISNULL) // empty cell has a chance to grow new tree

                        forest_n[i][j] = (ran1()                else // with fire, the tree is dead

                        forest_n[i][j] = ISNULL;

            }

        }

        // update

        for(i=0;i

            for(j=0;j

                forest[i][j]=forest_n[i][j];

        // store data

        for(i=0;i

            for(j=0;j

                fprintf(pFile,"%2d ",forest[i][j]);

            fprintf(pFile,"\n");

        }

    }

    fclose(pFile);

    return 0;

}

This post contains some useful links for c/c++ learning~~

1. Understanding the “extern” key word in C

http://www.geeksforgeeks.org/archives/840

(more to be added)

 

—-use profiler

[1]: g++ -o sortrand -pg main.cpp
[2]: ./sortrand
[3]: gprof sortrand.exe
—-显示库中的函数： ar t name.a

—-useful settings of output format in C++
cout.setf(ios_base::scientific,ios_base::floatfield);
cout.setf(ios_base::left,ios_base::adjustfield);
cout.width(0);恢复缺省设置

孤立子是数学和物理上最美丽的存在之一(待续)