题目：http://www.lydsy.com/JudgeOnline/problem.php?id=1065
这个图就是一个环套了一些树。
我们可以发现每个点对答案的贡献是(Ci k^d )/(1-k^len) len表示环的长度，然后这样的话我们某一个点修改一定是改到1下面去，那么考虑环上修改的点会让环长度改变，那么枚举这些点中离1最远的那个，然后断开环，弄成一棵内向树，然后在内向树上DP，状态dp(a,b,c)表示子树a，深度b，修改了c次的最大值，然后转移就分b=1&&dep(a)!=1（修改了a）和不修改a讨论即可，这样做的复杂度就是O(n^5)，常数比较小，所以可以AC。
代码：
#include <cstdio>
#include <algorithm>
#include <cstring>
using namespace std ;
#define REP( i , l , r ) for ( int i = l ; i <= r ; ++ i )
#define rep( i , x ) for ( int i = 0 ; i ++ < x ; )
#define DOWN( i , r , l ) for ( int i = r ; i >= l ; -- i )
const int maxn = 65 ;
inline void update( double &x , double y ) {
if ( y > x ) x = y ;
}
double k , c[ maxn ] , mul[ maxn ] , ans ;
int n , next[ maxn ] , m ;
bool flag[ maxn ] ;
struct edge {
edge *next ;
int t ;
} E[ maxn ] ;
edge *pt , *head[ maxn ] ;
inline void Init_edge( ) {
memset( head , 0 , sizeof( head ) ) , pt = E ;
}
inline void addedge( int s , int t ) {
pt -> t = t , pt -> next = head[ s ] ; head[ s ] = pt ++ ;
}
#define travel( x ) for ( edge *p = head[ x ] ; p ; p = p -> next )
double dp[ maxn ][ maxn ][ maxn ] ;
int sz[ maxn ] , ht[ maxn ] ;
void geth( int now ) {
sz[ now ] = 1 ;
travel( now ) {
ht[ p -> t ] = ht[ now ] + 1 ;
geth( p -> t ) ;
sz[ now ] += sz[ p -> t ] ;
}
}
const double inf = double( 100000000 ) * double( 100000000 ) ;
double g[ maxn ][ maxn ] ;
void dfs( int now ) {
REP( i , 0 , n ) REP( j , 0 , n ) dp[ now ][ i ][ j ] = - inf ;
if ( ! head[ now ] ) {
REP( i , 2 , ht[ now ] ) dp[ now ][ i ][ 0 ] = c[ now ] * mul[ i ] ;
if ( ht[ now ] != 1 ) dp[ now ][ 1 ][ 1 ] = c[ now ] * k ; else dp[ now ][ 1 ][ 0 ] = c[ now ] * k ;
return ;
}
travel( now ) dfs( p -> t ) ;
REP( dep , 0 , ht[ now ] ) {
if ( ! dep && now != 1 ) continue ;
int cnt = 0 ;
double cost ;
REP( i , 0 , sz[ now ] ) g[ 0 ][ i ] = - inf ; g[ 0 ][ 0 ] = 0 ;
travel( now ) {
++ cnt ;
REP( i , 0 , sz[ now ] ) g[ cnt ][ i ] = - inf ;
REP( i , 0 , sz[ p -> t ] ) {
cost = max( dp[ p -> t ][ dep + 1 ][ i ] , dp[ p -> t ][ 1 ][ i ] ) ;
REP( j , i , sz[ now ] ) {
update( g[ cnt ][ j ] , g[ cnt - 1 ][ j - i ] + cost ) ;
}
}
}
REP( tmp , 0 , sz[ now ] ) {
if ( dep == 1 && ht[ now ] != 1 ) {
if ( tmp ) dp[ now ][ dep ][ tmp ] = g[ cnt ][ tmp - 1 ] + c[ now ] * mul[ dep ] ;
} else {
dp[ now ][ dep ][ tmp ] = g[ cnt ][ tmp ] + c[ now ] * mul[ dep ] ;
}
}
}
}
inline double solve( ) {
int len = 0 ; for ( int i = 1 ; ; i = next[ i ] ) {
++ len ; if ( next[ i ] == 1 ) break ;
}
Init_edge( ) ;
REP( i , 2 , n ) addedge( next[ i ] , i ) ;
ht[ 1 ] = 0 ; geth( 1 ) ;
dfs( 1 ) ;
double tmp = - inf ;
REP( i , 0 , m ) update( tmp , dp[ 1 ][ 0 ][ i ] ) ;
return tmp / ( 1.0 - mul[ len ] ) ;
}
int main( ) {
scanf( "%d%d%lf" , &n , &m , &k ) ;
rep( i , n ) scanf( "%d" , next + i ) ; rep( i , n ) scanf( "%lf" , c + i ) ;
mul[ 0 ] = 1 ; rep( i , n ) mul[ i ] = mul[ i - 1 ] * k ;
ans = solve( ) ;
memset( flag , false , sizeof( flag ) ) ;
for ( int i = 1 ; ; i = next[ i ] ) {
flag[ i ] = true ;
if ( next[ i ] == 1 ) break ;
}
if ( m ) {
-- m ;
rep( i , n ) if ( flag[ i ] && i != 1 && next[ i ] != 1 ) {
int t = next[ i ] ; next[ i ] = 1 ;
update( ans , solve( ) ) ;
next[ i ] = t ;
}
}
printf( "%.2f\n" , ans ) ;
return 0 ;
}
